import numpy as np
import pandas as pd
import pandasql as ps
import math
import itertools 
import calendar
import time
from datetime import timedelta
import datetime

# matplotlib libraries
import matplotlib.pyplot as plt
from matplotlib.ticker import FuncFormatter
import matplotlib.cm as cm
import matplotlib.colors
import seaborn as sns

# date libraries
from dateutil import parser
from datetime import datetime, timedelta, date
import holidays

# prophet library
from prophet import Prophet
from prophet.diagnostics import performance_metrics
from prophet.plot import plot_cross_validation_metric
from prophet.diagnostics import cross_validation

# SKLearn libraries
from sklearn.metrics import mean_squared_error
from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import StandardScaler

import os
import warnings

# Suppress FutureWarnings specifically
warnings.filterwarnings("ignore", category=FutureWarning)

# Set Pandas option to display the entire content in a column
pd.set_option('display.max_columns', None)

# Define the width for metrics display alignment
label_width = 10
value_width = 15

def calculate_metrics(df):
    '''
    Function to calculate RMSE, MAPE, and Accuracy
    
    input: DataFrame with columns 'y' and 'yhat'
    output: RMSE, MAPE, and Accuracy values
    '''
    if 'y' not in df.columns or 'yhat' not in df.columns:
        raise KeyError("The DataFrame must contain 'y' and 'yhat' columns.")

    # Calculate RMSE - Root Mean Square Error
    rmse = np.sqrt(mean_squared_error(df['y'], df['yhat']))
    
    # Calculate MAE - Mean Absolute Error
    mae = np.mean(abs(df['y'] - df['yhat']))
    
    # Calculate MAPE - Mean Absolute Percentage Error
    mape = np.mean(abs((df['y'] - df['yhat']) / df['y'])) * 100
    
    # Calculate Accuracy
    accuracy = 100 - mape
    
    return rmse, mae, mape, accuracy

# Credit to The Almost Astrophysicist for this function
def missing_data(input_data):
    '''
    This function returns dataframe with information about the percentage of nulls in each column and the column data type.
    
    input: pandas df
    output: pandas df
    '''
    
    total = input_data.isnull().sum()
    percent = (input_data.isnull().sum()/input_data.isnull().count()*100)
    table = pd.concat([total, percent], axis = 1, keys = ['Total', 'Percent'])
    types = []
    for col in input_data.columns: 
        dtype = str(input_data[col].dtype)
        types.append(dtype)
    table["Types"] = types
    return(pd.DataFrame(table))

# Credit to Rob Mulla for this function
from pandas.api.types import CategoricalDtype

cat_type = CategoricalDtype(categories=['Monday', 'Tuesday', 'Wednesday',
                                        'Thursday', 'Friday', 'Saturday', 'Sunday'],
                            ordered=True)

def create_time_features(df, label=None):
    """
    Creates time-series features from datetime index and includes all original columns.
    Adds features for the first weekend of the month and the first weekend after the 15th.
    """
    print(f'The label is: {label}')
    df = df.copy()
    df['date'] = df.index
    df['hour'] = df['date'].dt.hour
    df['dayofweek'] = df['date'].dt.dayofweek
    df['weekday'] = df['date'].dt.day_name()
    df['weekday'] = df['weekday'].astype(cat_type)
    df['quarter'] = df['date'].dt.quarter
    df['month'] = df['date'].dt.month
    df['month_name'] = df['date'].dt.month_name()  # Add the name of the month
    df['year'] = df['date'].dt.year
    df['dayofyear'] = df['date'].dt.dayofyear
    df['dayofmonth'] = df['date'].dt.day
    df['weekofyear'] = df.index.isocalendar().week
    df['date_offset'] = (df.date.dt.month * 100 + df.date.dt.day - 320) % 1300

    """
    US seasons are not applicable in Ecuador
    df['season'] = pd.cut(df['date_offset'], [0, 300, 602, 900, 1300],
                          labels=['Spring', 'Summer', 'Fall', 'Winter'])
    """
    
    # Define wet and dry seasons based on Ecuador's climate
    df['season'] = df['month'].apply(lambda x: 'Wet' if x in [11, 12, 1, 2, 3, 4, 5] else 'Dry')
    
    # Create 'First Weekend of Month' feature
    df['first_weekend_of_month'] = 0
    for year, month in df[['year', 'month']].drop_duplicates().itertuples(index=False):
        month_df = df[(df['year'] == year) & (df['month'] == month)]
        if not month_df.empty:
            first_day = month_df.index.min()
            first_saturday = first_day + pd.DateOffset(days=(5 - first_day.dayofweek) % 7)
            first_sunday = first_day + pd.DateOffset(days=(6 - first_day.dayofweek) % 7)
            
            # Ensure the calculated dates exist within the month
            if first_saturday in df.index:
                df.loc[first_saturday, 'first_weekend_of_month'] = 1
            if first_sunday in df.index:
                df.loc[first_sunday, 'first_weekend_of_month'] = 1

    # Create 'First Weekend After 15th' feature
    df['first_weekend_after_15th'] = 0
    for year, month in df[['year', 'month']].drop_duplicates().itertuples(index=False):
        month_df = df[(df['year'] == year) & (df['month'] == month)]
        if not month_df.empty:
            fifteenth_day = month_df.loc[month_df['dayofmonth'] == 15].index[0]
            first_saturday_after_15th = fifteenth_day + pd.DateOffset(days=(5 - fifteenth_day.dayofweek) % 7)
            first_sunday_after_15th = fifteenth_day + pd.DateOffset(days=(6 - fifteenth_day.dayofweek) % 7)
            
            # Ensure the calculated dates exist within the month
            if first_saturday_after_15th in df.index:
                df.loc[first_saturday_after_15th, 'first_weekend_after_15th'] = 1
            if first_sunday_after_15th in df.index:
                df.loc[first_sunday_after_15th, 'first_weekend_after_15th'] = 1

    # Create 'Week Before Christmas' feature (18th to 24th December)
    #df['is_week_before_christmas'] = (df['month'] == 12) & (df['dayofmonth'].between(18, 24))
    
    # Return the DataFrame with all original columns plus the new features
    return df

def plot_train_set(test_start_date, df_pred_plt, model_name):
    # Create an explicit copy of the DataFrame to avoid warnings
    df_pred_plt = df_pred_plt.copy()

    # Convert the 'ds' column to datetime
    df_pred_plt['ds'] = pd.to_datetime(df_pred_plt['ds'])

    # Split data based on the test_start_date
    before_test = df_pred_plt[df_pred_plt['ds'] < test_start_date]
    on_or_after_test = df_pred_plt[df_pred_plt['ds'] >= test_start_date]

    # Plot the time series
    plt.figure(figsize=(14, 8))

    # Plot the actual sales
    plt.plot(df_pred_plt['ds'], df_pred_plt['y'], alpha=1, label='Actual Grocery Sales')

    # Add main title and subtitle
    plt.suptitle('Total Daily Sales of "Grocery I" for All Stores Combined', fontsize=16)
    plt.title(f'Predicted Sales vs Actual Sales (Model: {model_name}). 15 Aug 2015 - 15 Aug 2017', fontsize=12)

    plt.xlabel('Date')
    plt.ylabel('Total Daily Sales ($)')
    plt.grid(True)

    # Plot predicted sales before test_start_date in green
    plt.plot(before_test['ds'], before_test['yhat'], color='green', alpha=1, label='Predicted Sales (Training Set)')

    # Plot predicted sales on or after test_start_date in red
    plt.plot(on_or_after_test['ds'], on_or_after_test['yhat'], color='red', alpha=1, label='Predicted Sales (Test Set)')

    # Add legend
    plt.legend()
    plt.show()

def plot_test_set(model_name, df_pred_plt, df_test=None):
    # Define the new start date for the plot
    nov15_start_date = '2016-11-15'

    # Filter the DataFrame for dates on or after the new start date
    filtered_df = df_pred_plt[df_pred_plt['ds'] >= nov15_start_date]

    # Split data based on the test_start_date
    before_test = filtered_df[filtered_df['ds'] < test_start_date]

    # Check if df_test is provided
    if df_test is None:
        on_or_after_test = filtered_df[filtered_df['ds'] >= test_start_date]
    else:
        on_or_after_test = df_test

    # Plot the time series
    plt.figure(figsize=(14, 8))

    # Plot the actual sales
    plt.plot(filtered_df['ds'], filtered_df['y'], alpha=1, label='Actual Grocery Sales')

    # Add main title and subtitle
    plt.suptitle('Total Daily Sales of "Grocery I" for All Stores Combined', fontsize=16)
    plt.title(f'Predicted Sales vs Actual Sales (Model: {model_name}). 15 Nov 2016 - 15 Aug 2017', fontsize=12)

    plt.xlabel('Date', fontsize=14)
    plt.ylabel('Total Daily Sales ($)', fontsize=14)
    plt.grid(True)

    # Plot predicted sales before test_start_date in green
    plt.plot(before_test['ds'], before_test['yhat'], color='green', alpha=1, label='Predicted Sales (Training Set)')
    
    # Plot confidence interval for the test set
    plt.fill_between(before_test['ds'], before_test['yhat_lower'], before_test['yhat_upper'],
                     color='green', alpha=0.25, label='80% Confidence Interval (Training Set)')


    # Plot predicted sales on or after test_start_date in red
    plt.plot(on_or_after_test['ds'], on_or_after_test['yhat'], color='red', alpha=1, label='Predicted Sales (Test Set)')
    
    # Plot confidence interval for the test set
    plt.fill_between(on_or_after_test['ds'], on_or_after_test['yhat_lower'], on_or_after_test['yhat_upper'],
                     color='red', alpha=0.25, label='80% Confidence Interval (Test Set)')

    # Add legend
    plt.legend(fontsize=12)
    plt.show()

# Read the training dataset
df = pd.read_csv("./data/train.csv")

# read in the stores dataset
df_stores = pd.read_csv("./data/stores.csv")

# read in the stores dataset
df_holidays = pd.read_csv("./data/holidays_events.csv")

# read in the stores dataset
df_oil = pd.read_csv("./data/oil.csv")

df.head()

df.shape

(3000888, 6)

df['date'].nunique()

1684

df_stores.head()

df_holidays.head()

df_oil.head()

missing_data(df)

# data shape
print(df.shape)
# rows and columns.

(3000888, 6)

# summary of dataset
print(df.info())

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 3000888 entries, 0 to 3000887
Data columns (total 6 columns):
 #   Column       Dtype  
---  ------       -----  
 0   id           int64  
 1   date         object 
 2   store_nbr    int64  
 3   family       object 
 4   sales        float64
 5   onpromotion  int64  
dtypes: float64(1), int64(3), object(2)
memory usage: 137.4+ MB
None

df['store_nbr'].nunique(),

(54,)

min_date = df['date'].min()
unique_store_nbr_at_min_date = df[df['date'] == min_date]['store_nbr'].nunique()
print(unique_store_nbr_at_min_date)

54

df['family'].nunique()

33

min_date = df['date'].min()
unique_store_nbr_at_min_date = df[df['date'] == min_date]['family'].nunique()
print(unique_store_nbr_at_min_date)

33

# Convert the 'date' column to datetime format
df['date'] = pd.to_datetime(df['date'])

# Check if there are any NaT values, which indicate failed conversions
invalid_dates = df[df['date'].isna()]
invalid_dates.head()

print(df.info())

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 3000888 entries, 0 to 3000887
Data columns (total 6 columns):
 #   Column       Dtype         
---  ------       -----         
 0   id           int64         
 1   date         datetime64[ns]
 2   store_nbr    int64         
 3   family       object        
 4   sales        float64       
 5   onpromotion  int64         
dtypes: datetime64[ns](1), float64(1), int64(3), object(1)
memory usage: 137.4+ MB
None

# Join df with df_stores on 'store_nbr'
df = df.merge(df_stores, on='store_nbr', how='left')

# Check for missing data
missing_data(df_oil)

# Convert the date column into a datetime format
df_oil['date'] = pd.to_datetime(df_oil['date'], dayfirst=True, errors='coerce')

# Check if there are any NaT (Not a Time) values, which missing or null values in datetime columns.
invalid_dates = df_oil[df_oil['date'].isna()]
invalid_dates.head()

# Interpolate the missing values using linear interpolation
df_oil['dcoilwtico'] = df_oil['dcoilwtico'].interpolate(method='linear')

# Fill any remaining NaN values at the start or end using forward/backward fill
df_oil['dcoilwtico'] = df_oil['dcoilwtico'].ffill()
df_oil['dcoilwtico'] = df_oil['dcoilwtico'].bfill()

# Check for missing data
missing_data(df_oil)

# Join df with df_oil on 'date'
df = df.merge(df_oil, on='date', how='left')

# View the header
df.head()

df.info()

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 3000888 entries, 0 to 3000887
Data columns (total 11 columns):
 #   Column       Dtype         
---  ------       -----         
 0   id           int64         
 1   date         datetime64[ns]
 2   store_nbr    int64         
 3   family       object        
 4   sales        float64       
 5   onpromotion  int64         
 6   city         object        
 7   state        object        
 8   type         object        
 9   cluster      int64         
 10  dcoilwtico   float64       
dtypes: datetime64[ns](1), float64(2), int64(4), object(4)
memory usage: 251.8+ MB

# Group by 'date' and sum 'sales' column
daily_total_sales = df.groupby('date')['sales'].sum()

# Plot the time series
plt.figure(figsize=(14, 8))
daily_total_sales.plot()
plt.title('Daily Total Sales for All Stores Combined')
plt.xlabel('Date')
plt.ylabel('Total Sales ($ Millions)')
plt.grid(True)
plt.show()

# Group by 'date' and 'family' and sum the 'sales' column
# Credit to The Almost Astrophysicist for this code.
daily_total_sales_by_product = df.groupby(['date', 'family'])['sales'].sum().unstack()

# Make a monthly df and repeat the same plots
df_monthly = df.copy()

# Create a new column 'year_month' for aggregation
df_monthly['year_month'] = df_monthly['date'].dt.to_period('M')

# Group by 'year_month', 'store_nbr', and 'family' and aggregate 'sales' and 'onpromotion'
df_monthly = df_monthly.groupby(['year_month', 'store_nbr', 'family'], as_index=False).agg({
    'sales': 'sum',
    'onpromotion': 'sum'
})

# Set the 'year_month' column as a datetime object again:
df_monthly['year_month'] = df_monthly['year_month'].dt.to_timestamp()

df_monthly.tail()

# Group by 'date' and sum 'sales' column
daily_total_sales = df_monthly.groupby('year_month')['sales'].sum()

# Plot the time series
plt.figure(figsize=(14, 8))
daily_total_sales.plot()
plt.title('Monthly Total Sales for All Stores Combined')
plt.xlabel('Date')
plt.ylabel('Total Sales ($ Millions)')
plt.grid(True)
plt.show()

# Group by 'date' and 'family' and sum the 'sales' column
# Credit to The Almost Astrophysicist for this code.
monthly_total_sales = df_monthly.groupby(['year_month', 'family'])['sales'].sum().unstack()

# Plot the daily sale history for each product category
plt.figure(figsize=(14, 8))
for family in monthly_total_sales.columns:
    monthly_total_sales[family].plot(label=family)

plt.title('Combined Stores - Daily Total Sales by Product Category')
plt.xlabel('Date')
plt.ylabel('Total Sales $')
plt.legend(loc='upper center', bbox_to_anchor=(0.5, -0.15), shadow=True, ncol=2)
plt.grid(True)
plt.show()

daily_total_sales_by_product

daily_total_sales_by_product.shape

(1684, 33)

# Create a figure and axes for subplots
fig, axs = plt.subplots(len(daily_total_sales_by_product.columns), 1, figsize=(12, 6 * len(daily_total_sales_by_product.columns)))

# Iterate through each product and plot its daily total sales
for i, product in enumerate(daily_total_sales_by_product.columns):
    axs[i].plot(daily_total_sales_by_product[product])
    axs[i].set_title(f'Combined Stores - Daily Total Sales of {product}')
    axs[i].set_xlabel('Date')
    axs[i].set_ylabel('Total Sales $')
    axs[i].grid(True)

# Adjust layout to prevent overlap
plt.tight_layout()

# Show the plots
plt.show()

df_grocery = df[df['family'] == 'GROCERY I']
dfg = df_grocery.groupby('date').agg({
    'sales': 'sum', # Total sales
    'onpromotion': 'mean'  # The mean value of promotion across all stores
}).reset_index()

# Convert the 'date' column to datetime format
dfg['date'] = pd.to_datetime(dfg['date'])

# Rejoin oil onto the dataset
dfg = dfg.merge(df_oil, on='date', how='left')

missing_data(dfg)

# Fill missing oil prices by interpolation
dfg['dcoilwtico'] = dfg['dcoilwtico'].interpolate(method='linear')
dfg['dcoilwtico'] = dfg['dcoilwtico'].ffill()
dfg['dcoilwtico'] = dfg['dcoilwtico'].bfill()

# Check for missing data
missing_data(dfg)

# Ensure 'date' is the index
dfg.set_index('date', inplace=True)

# Plot the time series
plt.figure(figsize=(14, 8))
dfg['sales'].plot()
plt.title('Daily Total Sales of "Grocery I" for All Stores Combined')
plt.xlabel('Date')
plt.ylabel('Total Sales ($US)')
plt.grid(True)
plt.show()

print(dfg.index.min())  # Ensure the minimum date is 2013-01-01

2013-01-01 00:00:00

# Plot the oil prices
plt.figure(figsize=(14, 8))
dfg['dcoilwtico'].plot()
plt.title('Daily Oil Prices per Barrel')
plt.xlabel('Date')
plt.ylabel('Price per Barrel ($US)')
plt.grid(True)
plt.show()

# Setup the figure and axes
fig, ax = plt.subplots(figsize=(14, 8))

# Main Title

# Drop the rows where 'dcoilwtico' or 'sales' have NaN values
df_plot = dfg.dropna(subset=['dcoilwtico', 'sales'])
df_plot = df_plot[~(np.isinf(df_plot['dcoilwtico']) | np.isinf(df_plot['sales']))]

# Create a scatter plot for 'dcoilwtico' vs 'sales'
sns.scatterplot(x=df_plot['dcoilwtico'], y=df_plot['sales'], ax=ax)

# Create a line of best fit for 'dcoilwtico' vs 'sales'
m, b = np.polyfit(df_plot['dcoilwtico'], df_plot['sales'], 1)

ax.plot(df_plot['dcoilwtico'], m*df_plot['dcoilwtico'] + b, color='red')
ax.set_title("Daily Oil Prices per Barrel vs Sales", fontsize=18)
ax.set_xlabel("Daily Oil Prices per Barrel (US$)", fontsize=15)
ax.set_ylabel("Sales", fontsize=15)
ax.tick_params(axis='both', labelsize=12)

# Display the plot
plt.tight_layout()
plt.show()

# Plot the average number of promotions of groceries per store 
plt.figure(figsize=(14, 8))
dfg['onpromotion'].plot()
plt.title('Daily Average Number of Sale Promotions per Store')
plt.xlabel('Date')
plt.ylabel('Average number of Promotions per Store')
plt.grid(True)
plt.show()

# Plot the oil prices
plt.figure(figsize=(14, 8))
np.log1p(dfg['onpromotion']).plot()
plt.title('Natural Log of Daily Average Number of Sale Promotions per Store')
plt.xlabel('Date')
plt.ylabel('Natural Log Average number of Promotions per Store')
plt.grid(True)
plt.show()

# Create a new feature of log(onpromotion)
dfg['log_onpromotion'] = np.log1p(dfg['onpromotion'])

# This operation is equivalent to
# np.log(dfg['onpromotion'] + 1)

# Setup the figure and axes
fig, axes = plt.subplots(nrows=2, ncols=1, figsize=(16, 16))

# Adjust spacing between the main title and subplots
fig.subplots_adjust(top=0.94)  # Adjust the top margin to create space for the main title

# Add the main title
fig.suptitle('Comparison of Cost and Log of Cost', fontsize=20)

# Plot 1 
# Drop the rows where 'onpromotion' or 'sales' have NaN values
df_plot1 = dfg.dropna(subset=['onpromotion', 'sales'])
df_plot1 = df_plot1[~(np.isinf(df_plot1['onpromotion']) | np.isinf(df_plot1['sales']))]

# Create a scatter plot for cost vs rating number
sns.scatterplot(x=df_plot1['onpromotion'], y=df_plot1['sales'], ax=axes[0])

# Create a line of best fit for cost vs rating number
m1, b1 = np.polyfit(df_plot1['onpromotion'], df_plot1['sales'], 1)
axes[0].plot(df_plot1['onpromotion'], m1*df_plot1['onpromotion'] + b1, color='red')

# Customize the first plot
axes[0].set_title("Average Number of Promotions per Store vs Sales", fontsize=18)
axes[0].set_xlabel("Average Number of Promotions per Store", fontsize=15)
axes[0].set_ylabel("Sales", fontsize=15)
axes[0].tick_params(axis='both', labelsize=12)

# Plot 2
# Drop the rows where 'log_onpromotion' or 'sales' have NaN values
df_plot2 = dfg.dropna(subset=['log_onpromotion', 'sales'])
df_plot2 = df_plot2[~(np.isinf(df_plot2['log_onpromotion']) | np.isinf(df_plot2['sales']))]

# Create a scatter plot for log_onpromotion vs sales
sns.scatterplot(x=df_plot2['log_onpromotion'], y=df_plot2['sales'], ax=axes[1])

# Create a line of best fit for log_onpromotion vs sales
m2, b2 = np.polyfit(df_plot2['log_onpromotion'], df_plot2['sales'], 1)
axes[1].plot(df_plot2['log_onpromotion'], m2*df_plot2['log_onpromotion'] + b2, color='red')

# Customize the second plot
axes[1].set_title("Log of Average Number of Promotions per Store vs Sales", fontsize=18)
axes[1].set_xlabel("Log of Average Number of Promotions per Store", fontsize=15)
axes[1].set_ylabel("Sales", fontsize=15)
axes[1].tick_params(axis='both', labelsize=12)

# Display the plots
plt.show()

dfg.info()

<class 'pandas.core.frame.DataFrame'>
DatetimeIndex: 1684 entries, 2013-01-01 to 2017-08-15
Data columns (total 4 columns):
 #   Column           Non-Null Count  Dtype  
---  ------           --------------  -----  
 0   sales            1684 non-null   float64
 1   onpromotion      1684 non-null   float64
 2   dcoilwtico       1684 non-null   float64
 3   log_onpromotion  1684 non-null   float64
dtypes: float64(4)
memory usage: 65.8 KB

dfg.describe()

dfg.shape

(1684, 4)

full_date_range = pd.date_range(start="2013-01-01", end="2017-08-15")
missing_dates = full_date_range.difference(dfg.index)
print("Missing dates:", missing_dates)

Missing dates: DatetimeIndex(['2013-12-25', '2014-12-25', '2015-12-25', '2016-12-25'], dtype='datetime64[ns]', freq=None)

plt.figure(figsize=(14, 8))
dfg['sales'].plot(
    kind='hist',
    bins=200,
    title='Distribution of Daily Sales',
    xlabel='Daily Grocery Sales ($US)',
    ylabel='Number of Days'
)
plt.show()

# Determine the date of maximum sales in the groceries and the quantity
dfg[dfg['sales'] == dfg['sales'].max()]

# Call the create_time_features function to create time features for the data frame
dfgt = create_time_features(dfg)

The label is: None

# Make a copy of the time featured dataset to be used for the prophet model.
dfg1 = dfgt.copy()

dfg1.head(6)

# Get the month names
month_names = [calendar.month_name[i] for i in range(1, 13)]

# Create the plot
fig, ax = plt.subplots(figsize=(14, 8))
sns.boxplot(data=dfg1.dropna(),
            x='month',
            y='sales',
            hue='weekday',
            ax=ax,
            linewidth=1)

# Set the x-axis labels to the month names
ax.set_xticklabels(month_names, rotation=45, ha='right')

# Set the title, x-label, y-label, and legend
ax.set_title('Box Plot Grocery Sales by Month and Day of Week')
ax.set_xlabel('Month')
ax.set_ylabel('Daily Grocery Sales ($)')
ax.legend(bbox_to_anchor=(1, 1))

# Show the plot
plt.show()

# Select a data frame for the first 6 months of 2014
df_2014_1_3 = dfg1.loc[('2014-01-01' <= dfg1.index) & (dfg1.index <= '2014-06-30')]

# Plot the sales for 01 January to 30 June 2014
plt.figure(figsize=(14, 8))
df_2014_1_3['sales'].plot()
plt.title('Daily Total Sales of "Grocery I" for All Stores Combined in January to June 2014')
plt.xlabel('Date')
plt.ylabel('Total Sales ($)')
plt.grid(True)
plt.show()

# Define the test start date
test_start_date = "2017-01-01"
test_start_date = datetime.strptime(test_start_date, "%Y-%m-%d")

# Convert test_start_date to a Timestamp object
test_start_date_ts = pd.Timestamp(test_start_date)

# Make a copy of df grocery (dfg1) to df Prophet (dfp)
dfp = dfg1.copy()

dfp.head()

# Step 1: Drop the existing 'date' column if it exists
if 'date' in dfp.columns:
    dfp = dfp.drop(columns=['date'])

# Step 2: Reset the index to create the 'ds' column
dfp = dfp.reset_index()

# Step 3: Rename the columns: 'sales' to 'y' and 'index' (former date index) to 'ds'
dfp = dfp.rename(columns={'sales': 'y', 'date': 'ds'})

# Step 4: Keep only the 'ds' and 'y' columns
dfp = dfp[['ds', 'y']]

# Step 5: Ensure 'y' is numeric
dfp['y'] = pd.to_numeric(dfp['y'], errors='coerce')

# Display the Prophet Ready dataset information
dfp.info()

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 1684 entries, 0 to 1683
Data columns (total 2 columns):
 #   Column  Non-Null Count  Dtype         
---  ------  --------------  -----         
 0   ds      1684 non-null   datetime64[ns]
 1   y       1684 non-null   float64       
dtypes: datetime64[ns](1), float64(1)
memory usage: 26.4 KB

df_train = dfp[(dfp['ds'] < test_start_date) ]
#df_test = df_pred[df_pred['ds'] >= test_start_date]

# View the train tail
df_train.tail()

# Create an instance of the model and add the holiday data to it
pr1 = Prophet()

# Train the model
pr1.fit(df_train)

13:11:11 - cmdstanpy - INFO - Chain [1] start processing
13:11:12 - cmdstanpy - INFO - Chain [1] done processing

<prophet.forecaster.Prophet at 0x1c097d88e80>

dfg.tail()

# Find the last date in the data frame
last_date = dfg1['date'].max()

# Convert last_date to np.datetime64 with 'D' (day) unit
last_date_np = np.datetime64(last_date, 'D')

# Convert test_start_date_ts to np.datetime64 with 'D' (day) unit
test_start_date_np = np.datetime64(test_start_date_ts, 'D')

# Calculate the number of days between last_date and test_start_date + 1 to include the first day
num_days = np.timedelta64(last_date_np - test_start_date_np, 'D').astype(int) +1

num_days

227

# Make a future data frame for the number of specified days. 
# This is seven and a half months
future = pr1.make_future_dataframe(periods=num_days)

# Predict sales for the training dataset and the specified future number of days
df_pred = pr1.predict(future)

df_pred.head()

plt.figure(figsize=(14, 8))
fig1 = pr1.plot(df_pred)
plt.show()

<Figure size 1400x800 with 0 Axes>

# Plot the time series
plt.figure(figsize=(14, 8))
dfg['sales'].plot()
plt.title('Daily Total Sales of "Grocery I" for All Stores Combined')
plt.xlabel('Date')
plt.ylabel('Total Sales ($)')
plt.grid(True)
plt.show()

plt.figure(figsize=(14, 8))
pr1.plot_components(df_pred)
plt.show()

<Figure size 1400x800 with 0 Axes>

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# We are mainly interested in the predictions including and beyond the cutoff date of 2017-01-01
df_test1 = df_pred[df_pred['ds'] >= test_start_date]

df_test1 = df_test1.merge(dfp)

df_test1.head()

# Merge predictions with actual data over the entire dataset
df_pred1_long = df_pred.merge(dfp)

# Create a plot dataset showing the last two years of the predictions.
# This comprises both the train and test set.
df_pred1_plt = df_pred1_long[df_pred1_long['ds'] >= '2015-08-15']

model_name = 'Pr1'
plot_train_set(test_start_date, df_pred1_plt, model_name)

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# Filter for training data only, before the cutoff date of 2017-01-01
df_train1 = df_pred[df_pred['ds'] < test_start_date]

# Merge the actual sales values to the predicted test set
df_train1 = df_train1.merge(dfp)

# Call the function to calculate the performance metrics
rmse_trn1, mae_trn1, mape_trn1, accuracy_trn1 = calculate_metrics(df_train1)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_trn1, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_trn1, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_trn1, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_trn1, 2)}%')

RMSE:      31455.24       
MAE:       20612.7        
MAPE:      39.46%
Accuracy:  60.54%

plot_test_set(model_name, df_pred1_plt)

# Call the function to calculate the performance metrics
rmse_pr1, mae_pr1, mape_pr1, accuracy_pr1 = calculate_metrics(df_test1)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_pr1, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_pr1, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_pr1, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_pr1, 2)}%')

RMSE:      41012.74       
MAE:       24737.98       
MAPE:      81.51%
Accuracy:  18.49%

# Create a column of the daily error
df_test1['Error'] = df_test1['yhat'] - df_test1['y']
df_test1['Abs_Error'] = np.abs(df_test1['yhat'] - df_test1['y'])

df_temp = df_test1[['ds','y','yhat','Error','Abs_Error']]

df_temp.nlargest(5,'Abs_Error')

# Add the holiday data
holiday = pd.DataFrame([])

# Fetch Ecuadorian holidays from 2013 to 2017 inclusive
ecuador_holidays = holidays.Ecuador(years=[2013, 2014, 2015, 2016, 2017])

# Convert the holiday dates to a DataFrame
# It was determined that there is a strong sales effect one day before the holiday and two days after.
# These final numbers were confirmed, after evaluating many iterations.
ecuador_holidays_df = pd.DataFrame({
    'holiday': list(ecuador_holidays.values()),  # Name of the holiday
    'ds': pd.to_datetime(list(ecuador_holidays.keys())),  # Date of the holiday
    'lower_window': -1,  # Lower window
    'upper_window': 2,   # Upper window
})

# Ensure the lower and upper window for "Año Nuevo" (New Year) are set to 0
# As this holiday has virtually 0 sales, the effects in the lower and upper window were different to other holidays.
ecuador_holidays_df.loc[ecuador_holidays_df['holiday'] == 'Año Nuevo', ['lower_window', 'upper_window']] = [0, 0]
# Christmas day was briefly implemented into the dataset, but it resulted difficult to model, like New Year's so it was removed.
#ecuador_holidays_df.loc[ecuador_holidays_df['holiday'] == 'Día de Navidad', ['lower_window', 'upper_window']] = [0, 0]

# Display the holidays for the first two years
ecuador_holidays_df.sort_values('ds').head(25)

# Display the number of holidays in Ecuador during the time period.
ecuador_holidays_df.shape

(58, 4)

# Create an instance of the model and add the holiday data to it
pr2 = Prophet(holidays=ecuador_holidays_df)#,
             #seasonality_mode='additive')

# Train the model
pr2.fit(df_train)

13:11:14 - cmdstanpy - INFO - Chain [1] start processing
13:11:15 - cmdstanpy - INFO - Chain [1] done processing

<prophet.forecaster.Prophet at 0x1c097d91360>

# The training data consists of 4 years of daily sales, as shown below:
df_train

# Make a future data frame for the number of specified days. 
# This is seven and a half months
future2 = pr2.make_future_dataframe(periods=num_days)

# Predict sales for the training dataset and the specified future number of days
df_pred2 = pr2.predict(future2)

df_pred2.head()

plt.figure(figsize=(14, 8))
fig2 = pr2.plot(df_pred2)
plt.show()

<Figure size 1400x800 with 0 Axes>

plt.figure(figsize=(14, 8))
pr2.plot_components(df_pred2)
plt.show()

<Figure size 1400x800 with 0 Axes>

# View the holidays in the training set
pr2.train_holiday_names

0                                     Año Nuevo
1                                      Carnaval
2                                 Viernes Santo
3                               Día del Trabajo
4                          Batalla de Pichincha
5                 Primer Grito de Independencia
6                    Independencia de Guayaquil
7                           Día de los Difuntos
8                       Independencia de Cuenca
9                                Día de Navidad
10                        Año Nuevo (observado)
11             Batalla de Pichincha (observado)
12    Primer Grito de Independencia (observado)
dtype: object

df_pred2[df_pred2["Año Nuevo"] != 0]["Año Nuevo"]

0      -201583.851542
364    -201583.851542
728    -201583.851542
1092   -201583.851542
1457   -201583.851542
Name: Año Nuevo, dtype: float64

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# For now, we are mainly interested in the predictions including and beyond the cutoff date of 2017-06-15
df_test2 = df_pred2[df_pred2['ds'] >= test_start_date]

# Merge the actual values to the predictions in the test set.
df_test2 = df_test2.merge(dfp)

# Merge predictions with actual data over the entire dataset
df_pred2_long = df_pred2.merge(dfp)

# Create a plot dataset showing the last two years of the predictions.
# This comprises both the train and test set.
df_pred2_plt = df_pred2_long[df_pred2_long['ds'] >= '2015-08-15']

df_pred2_long.head()

model_name = 'Pr2.1'
plot_train_set(test_start_date, df_pred2_plt, model_name)

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# Filter for training data only, before the cutoff date of 2017-01-01
df_train2 = df_pred2[df_pred2['ds'] < test_start_date]

# Merge the actual sales values to the predicted test set
df_train2 = df_train2.merge(dfp)

# Call the function to calculate the performance metrics
rmse_trn2, mae_trn2, mape_trn2, accuracy_trn2 = calculate_metrics(df_train2)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_trn2, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_trn2, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_trn2, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_trn2, 2)}%')

RMSE:      27616.41       
MAE:       18490.93       
MAPE:      13.53%
Accuracy:  86.47%

plot_test_set(model_name, df_pred2_plt)

# Call the function to calculate the performance metrics
rmse_pr2, mae_pr2, mape_pr2, accuracy_pr2 = calculate_metrics(df_test2)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_pr2, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_pr2, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_pr2, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_pr2, 2)}%')

RMSE:      33199.77       
MAE:       21876.56       
MAPE:      40.2%
Accuracy:  59.8%

# Create a column of the daily error
df_test2['Error'] = df_test2['yhat'] - df_test2['y']
df_test2['Abs_Error'] = np.abs(df_test2['yhat'] - df_test2['y'])

df_temp = df_test2[['ds','y','yhat','Error','Abs_Error']]

df_temp.nlargest(5,'Abs_Error')

# Create an additive holidays data frame the excludes Christmas and New Year's Day.
additive_holidays_df = ecuador_holidays_df[
    ~ecuador_holidays_df['holiday'].isin(['Día de Navidad', 'Año Nuevo'])
]

# Display the holidays for the first two years
ecuador_holidays_df.sort_values('ds').head(15)

additive_holidays_df.head(15)

# Create binary flags for multiplicative holidays
dfg1['Día de Navidad'] = (dfg1.index.month == 12) & (dfg1.index.day == 25)
dfg1['Año Nuevo'] = (dfg1.index.month == 1) & (dfg1.index.day == 1)

# Convert to integers
dfg1['Día de Navidad'] = dfg1['Día de Navidad'].astype(int)
dfg1['Año Nuevo'] = dfg1['Año Nuevo'].astype(int)

dfg1.head()

ecuador_holidays_df.shape, additive_holidays_df.shape

((58, 4), (48, 4))

# Make a copy of df grocery (dfg) to df predict (dfp)
dfp = dfg1.copy()

# Prepare the dataframe for modelling

# Step 1: Drop the existing 'date' column
# The index is the date. The date column should be dropped, to avoid duplication.
if 'date' in dfp.columns:
    dfp = dfp.drop(columns=['date'])

# Step 2: Reset the index to create the 'date' column
dfp = dfp.reset_index()

# Step 3: Rename the columns: 'sales' to 'y' and 'index' (former date index) to 'ds' as is required for Prophet
dfp = dfp.rename(columns={'sales': 'y', 'date': 'ds'})

# Step 4: Keep only the 'ds', 'y' and the chosen holiday regressor columns
dfp = dfp[['ds', 'y', 'Día de Navidad', 'Año Nuevo']]

# Step 5: Ensure 'y' is numeric
dfp['y'] = pd.to_numeric(dfp['y'], errors='coerce')

df_train = dfp[(dfp['ds'] < test_start_date) ]

df_train.head()

# Ensure additive_holidays_df is a standalone copy
additive_holidays_df = additive_holidays_df.copy()

# Create an instance of the model and add the holiday data to it
pr22 = Prophet(holidays=additive_holidays_df,
             seasonality_mode='additive')

# Add multiplicative regressors for Christmas and New Year
pr22.add_regressor('Día de Navidad', mode='multiplicative')
pr22.add_regressor('Año Nuevo', mode='multiplicative')

<prophet.forecaster.Prophet at 0x1c09ac0bc40>

# Train the model
pr22.fit(df_train)

13:11:17 - cmdstanpy - INFO - Chain [1] start processing
13:11:18 - cmdstanpy - INFO - Chain [1] done processing

<prophet.forecaster.Prophet at 0x1c09ac0bc40>

# Make a future data frame for the number of specified days. 
# This is seven and a half months
future22 = pr22.make_future_dataframe(periods=num_days)

# Merge the regressor columns onto the future dataframe
future22 = future22.merge(dfp[['ds', 'Día de Navidad', 'Año Nuevo']], on='ds', how='left')

# Predict sales for the next two months
df_pred22 = pr22.predict(future22)

df_pred22.head()

plt.figure(figsize=(14, 8))
fig22 = pr22.plot(df_pred22)
plt.show()

<Figure size 1400x800 with 0 Axes>

plt.figure(figsize=(14, 8))
pr22.plot_components(df_pred22)
plt.show()

<Figure size 1400x800 with 0 Axes>

# Show the multiplicative factor for the new year regressor.
df_pred22[df_pred22["Año Nuevo"] != 0]["Año Nuevo"]

0      -1.08038
364    -1.08038
728    -1.08038
1092   -1.08038
1457   -1.08038
Name: Año Nuevo, dtype: float64

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# Filter out the test set
df_test22 = df_pred22[df_pred22['ds'] >= test_start_date]

# Merge the actual values to the predictions in the test set.
df_test22 = df_test22.merge(dfp, on='ds', how='left')

# Merge predictions with actual data over the entire dataset
df_pred22_long = df_pred22.merge(dfp, on='ds', how='left')

# Create a plot dataset showing the last two years of the predictions.
# This comprises both the train and test set.
df_pred22_plt = df_pred22_long[df_pred22_long['ds'] >= '2015-08-15']

model_name = 'Pr2.2'
plot_train_set(test_start_date, df_pred22_plt, model_name)

print(f"'2013-01-01': {df_pred22['yhat'][df_pred2['ds'] == '2013-01-01']}")
print(f"'2014-01-01': {df_pred22['yhat'][df_pred2['ds'] == '2014-01-01']}")
print(f"'2015-01-01': {df_pred22['yhat'][df_pred2['ds'] == '2015-01-01']}")
print(f"'2016-01-01': {df_pred22['yhat'][df_pred2['ds'] == '2016-01-01']}")
print(f"'2017-01-01': {df_pred22['yhat'][df_pred2['ds'] == '2017-01-01']}")

'2013-01-01': 0    7696.824606
Name: yhat, dtype: float64
'2014-01-01': 364    6895.44005
Name: yhat, dtype: float64
'2015-01-01': 728   -12693.858961
Name: yhat, dtype: float64
'2016-01-01': 1092    1380.799456
Name: yhat, dtype: float64
'2017-01-01': 1457    83090.0744
Name: yhat, dtype: float64

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.

# Merge the actual sales values to the predicted test set
df_train22 = df_pred22[['ds', 'yhat']].merge(dfp, on='ds', how='left')

# Filter for training data only, before the cutoff date of 2017-01-01
df_train22 = df_train22[df_train22['ds'] < test_start_date]

# Call the function to calculate the performance metrics
rmse_trn22, mae_trn22, mape_trn22, accuracy_trn22 = calculate_metrics(df_train22)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_trn22, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_trn22, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_trn22, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_trn22, 2)}%')

RMSE:      28053.68       
MAE:       18671.96       
MAPE:      10.3%
Accuracy:  89.7%

plot_test_set(model_name, df_pred22_plt)

# Call the function to calculate the performance metrics
rmse_pr2, mae_pr2, mape_pr2, accuracy_pr2 = calculate_metrics(df_test22)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_pr2, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_pr2, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_pr2, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_pr2, 2)}%')

RMSE:      32534.5        
MAE:       21782.73       
MAPE:      25.66%
Accuracy:  74.34%

def model_new_years_day_prophet(df, df_test):
    # Filter data for January 1st only in the training set
    df_jan1 = df[(df['ds'].dt.month == 1) & (df['ds'].dt.day == 1)]

    # Create a new Prophet model
    model_jan1 = Prophet(yearly_seasonality=False, weekly_seasonality=False, daily_seasonality=False)
    
    # Fit the model on the January 1st data
    model_jan1.fit(df_jan1[['ds', 'y']])
    
    # Filter the test set for January 1st only
    df_test_jan1 = df_test[(df_test['ds'].dt.month == 1) & (df_test['ds'].dt.day == 1)]

    # Make predictions on the test set
    df_test_jan1 = model_jan1.predict(df_test_jan1[['ds']])

    # Insert the predictions back into the original test set
    df_test.loc[(df_test['ds'].dt.month == 1) & (df_test['ds'].dt.day == 1), 'yhat'] = df_test_jan1['yhat'].values
    
    # Plotting the train set and test set
    plt.figure(figsize=(14, 8))
    # Plot the training data
    plt.plot(df_jan1['ds'], df_jan1['y'], label='Train Set - Actual', marker='o')
    # Plot the test data and predictions
    plt.plot(df_test_jan1['ds'], df_test_jan1['yhat'], label='Test Set - Predicted', marker='x')
    if 'y' in df_test.columns:
        plt.plot(df_test_jan1['ds'], df_test.loc[df_test_jan1.index, 'y'], label='Test Set - Actual', marker='o')
    
    plt.title("Sales on January 1st - Train and Test Sets")
    plt.xlabel("Date")
    plt.ylabel("Sales")
    plt.legend()
    plt.grid(True)
    plt.show()
    
    return df_test

df_test2 = model_new_years_day_prophet(df_train, df_test2)

13:11:20 - cmdstanpy - INFO - Chain [1] start processing
13:11:20 - cmdstanpy - INFO - Chain [1] done processing

df_test2

model_name = 'Pr2.3'
plot_test_set(model_name, df_pred2_plt, df_test2)

# Call the function to calculate the performance metrics
rmse_pr2, mae_pr2, mape_pr2, accuracy_pr2 = calculate_metrics(df_test2)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_pr2, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_pr2, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_pr2, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_pr2, 2)}%')

RMSE:      31729.55       
MAE:       21238.83       
MAPE:      8.28%
Accuracy:  91.72%

# Create a column of the daily error
df_test2['Error'] = df_test2['yhat'] - df_test2['y']
df_test2['Abs_Error'] = np.abs(df_test2['yhat'] - df_test2['y'])

df_temp = df_test2[['ds','y','yhat','Error','Abs_Error']]

df_temp.nlargest(10,'Abs_Error')

# Select a data frame for the specified period in 2016
df_2016 = dfg1.loc[('2016-04-01' <= dfg1.index) & (dfg1.index <= '2016-05-05')]

# Plot the sales for April to May 2016
plt.figure(figsize=(14, 8))
plt.plot(df_2016.index, df_2016['sales'], label='Daily Sales')

# Add a red vertical line for the earthquake date
earthquake_date = '2016-04-16'
plt.axvline(pd.to_datetime(earthquake_date), color='red', linestyle='--', label='Manabi Earthquake (16 Apr 2016)')

# Add labels and title
plt.title('Daily Total Sales of "Grocery I" for All Stores Combined (April-May 2016)', fontsize=14)
plt.xlabel('Date', fontsize=12)
plt.xticks(rotation=45)
plt.ylabel('Total Sales ($)', fontsize=12)
plt.legend(fontsize=12)
plt.grid(True)

# Display the plot
plt.show()

max_sales_date = dfg1.loc[dfg1['sales'].idxmax(), 'date']
max_sales_value = dfg1['sales'].max()

print(f"The maximum sales value is: ${max_sales_value}")
print(f'It occurred on: {max_sales_date.strftime("%Y-%m-%d")}')

The maximum sales value is: $583222.48
It occurred on: 2016-04-18

# Define the specific dates
earthquake_1_day_after = '2016-04-17'
eq_2_days_after = '2016-04-18'
eq_3_days_after = '2016-04-19'
eq_4_days_after = '2016-04-20'
eq_5_days_after = '2016-04-21'
eq_6_days_after = '2016-04-22'

# Initialize the columns with 0s
dfg1['earthquake_1_day_after'] = 0
dfg1['eq_2_days_after'] = 0
dfg1['eq_3_days_after'] = 0
dfg1['eq_4_days_after'] = 0
dfg1['eq_5_days_after'] = 0
dfg1['eq_6_days_after'] = 0

# Set the flag for the specific dates with highly unusual activity
dfg1.loc[dfg1.index == earthquake_1_day_after, 'earthquake_1_day_after'] = 1
dfg1.loc[dfg1.index == eq_2_days_after, 'eq_2_days_after'] = 1
dfg1.loc[dfg1.index == eq_3_days_after, 'eq_3_days_after'] = 1
dfg1.loc[dfg1.index == eq_4_days_after, 'eq_4_days_after'] = 1
dfg1.loc[dfg1.index == eq_5_days_after, 'eq_5_days_after'] = 1
dfg1.loc[dfg1.index == eq_6_days_after, 'eq_6_days_after'] = 1
#dfg1.loc[dfg1.index == eq_7_days_after, 'eq_7_days_after'] = 1

dfg1.head()

# Make a copy of df grocery (dfg) to df predict (dfp)
dfp = dfg1.copy()

# Filter rows where either first_weekend_of_month or first_weekend_after_15th is 1
filtered_dfp = dfp[(dfp['first_weekend_of_month'] == 1) | (dfp['first_weekend_after_15th'] == 1)]

# Count the number of rows in the filtered DataFrame
num_rows = filtered_dfp.shape[0]  # This gives the number of rows

# Alternatively, you can use len() to count the rows
num_rows_alternative = len(filtered_dfp)

print(f"Number of days where weekend following payday regressors are present: {num_rows}")

Number of days where weekend following payday regressors are present: 222

# Step 1: Drop the existing 'date' column
# The index is the date. The date column should be dropped, to avoid duplication.
if 'date' in dfp.columns:
    dfp = dfp.drop(columns=['date'])

# Step 2: Reset the index to create the 'date' column
dfp = dfp.reset_index()

# Step 3: Rename the columns: 'sales' to 'y' and 'index' (former date index) to 'ds' as is required for Prophet
dfp = dfp.rename(columns={'sales': 'y', 'date': 'ds'})

# Step 4: Keep only the 'ds', 'y' and the chosen regressor columns
dfp = dfp[['ds', 'y', 'Año Nuevo', 'first_weekend_of_month', 'first_weekend_after_15th', 'dcoilwtico', 
           'onpromotion', 'log_onpromotion', 'earthquake_1_day_after', 'eq_2_days_after', 'eq_3_days_after', 
           'eq_4_days_after', 'eq_5_days_after', 'eq_6_days_after']]

# Step 5: Ensure 'y' is numeric
dfp['y'] = pd.to_numeric(dfp['y'], errors='coerce')

# Display the modified dataframe's info
dfp.info()

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 1684 entries, 0 to 1683
Data columns (total 14 columns):
 #   Column                    Non-Null Count  Dtype         
---  ------                    --------------  -----         
 0   ds                        1684 non-null   datetime64[ns]
 1   y                         1684 non-null   float64       
 2   Año Nuevo                 1684 non-null   int32         
 3   first_weekend_of_month    1684 non-null   int64         
 4   first_weekend_after_15th  1684 non-null   int64         
 5   dcoilwtico                1684 non-null   float64       
 6   onpromotion               1684 non-null   float64       
 7   log_onpromotion           1684 non-null   float64       
 8   earthquake_1_day_after    1684 non-null   int64         
 9   eq_2_days_after           1684 non-null   int64         
 10  eq_3_days_after           1684 non-null   int64         
 11  eq_4_days_after           1684 non-null   int64         
 12  eq_5_days_after           1684 non-null   int64         
 13  eq_6_days_after           1684 non-null   int64         
dtypes: datetime64[ns](1), float64(4), int32(1), int64(8)
memory usage: 177.7 KB

# Check for missing data
missing_data(dfp)

df_train3 = dfp[(dfp['ds'] < test_start_date) ]
#df_test = df_pred[df_pred['ds'] >= test_start_date]

df_train3.head()

# Create an instance of the model and add the holiday data to it
pr3 = Prophet(holidays=ecuador_holidays_df,
              seasonality_mode='additive')

# Add the regressors to the model
#pr3.add_regressor('Año Nuevo', mode='multiplicative')
pr3.add_regressor('first_weekend_of_month', mode='multiplicative')
pr3.add_regressor('first_weekend_after_15th', mode='multiplicative')
pr3.add_regressor('dcoilwtico')
#pr3.add_regressor('onpromotion')
#pr3.add_regressor('log_onpromotion')
pr3.add_regressor('earthquake_1_day_after')
pr3.add_regressor('eq_2_days_after')
pr3.add_regressor('eq_3_days_after')
pr3.add_regressor('eq_4_days_after')
pr3.add_regressor('eq_5_days_after')
pr3.add_regressor('eq_6_days_after')

<prophet.forecaster.Prophet at 0x1c09c14b6d0>

# Train the model
pr3.fit(df_train3)

13:11:22 - cmdstanpy - INFO - Chain [1] start processing
13:11:22 - cmdstanpy - INFO - Chain [1] done processing

<prophet.forecaster.Prophet at 0x1c09c14b6d0>

# Make a future data frame for the number of specified days. 
# This is seven and a half months
future3 = pr3.make_future_dataframe(periods=num_days)

# Merge the regressor columns onto the future dataframe
# Note: if they are commented out above, they will not be active even though they are included here.
future3 = future3.merge(dfp[['ds', 'Año Nuevo', 'first_weekend_of_month', 'first_weekend_after_15th', 'dcoilwtico', 
                             'onpromotion', 'log_onpromotion', 'earthquake_1_day_after', 'eq_2_days_after', 
                             'eq_3_days_after', 'eq_4_days_after', 'eq_5_days_after', 'eq_6_days_after']], 
                        on='ds', 
                        how='left')

# Predict sales for the next two months
df_pred3 = pr3.predict(future3)

plt.figure(figsize=(14, 8))
fig1 = pr3.plot(df_pred3)
plt.show()

<Figure size 1400x800 with 0 Axes>

plt.figure(figsize=(14, 8))
pr3.plot_components(df_pred3)
plt.show()

<Figure size 1400x800 with 0 Axes>

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# For now, we are only interested in the predictions including and beyond the cutoff date of 2017-01-01
df_test3 = df_pred3[df_pred3['ds'] >= test_start_date]

# Ensure df_test3 is a standalone copy
df_test3 = df_test3.copy()

# Convert 'ds' to datetime in both dataframes to ensure consistency
df_test3['ds'] = pd.to_datetime(df_test3['ds'])
dfp['ds'] = pd.to_datetime(dfp['ds'])

df_test3 = df_test3.merge(dfp[['ds', 'y']], on='ds', how='left')

df_test3.head()

# Merge predictions with actual data over the entire dataset
df_pred3_long = df_pred3.merge(dfp[['ds', 'y']], on='ds', how='left')

# Create a plot dataset showing the last two years of the predictions.
# This comprises both the train and test set.
df_pred3_plt = df_pred3_long[df_pred3_long['ds'] >= '2015-08-15']

df_pred3_long.tail()

model_name = 'Pr3.1'
plot_train_set(test_start_date, df_pred3_plt, model_name)

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.

# Merge the actual sales values to the predicted test set
df_train3 = df_pred3[['ds', 'yhat']].merge(dfp, on='ds', how='left')

# Filter for training data only, before the cutoff date of 2017-01-01
df_train3 = df_train3[df_train3['ds'] < test_start_date]

# Call the function to calculate the performance metrics
rmse_trn3, mae_trn3, mape_trn3, accuracy_trn3 = calculate_metrics(df_train3)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_trn3, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_trn3, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_trn3, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_trn3,2)}%')

RMSE:      22418.38       
MAE:       16197.24       
MAPE:      12.47%
Accuracy:  87.53%

df_test3 = model_new_years_day_prophet(df_train, df_test3)

13:11:25 - cmdstanpy - INFO - Chain [1] start processing
13:11:25 - cmdstanpy - INFO - Chain [1] done processing

plot_test_set(model_name, df_pred3_plt, df_test3)

# Call the function to calculate the performance metrics
rmse_pr3, mae_pr3, mape_pr3, accuracy_pr3 = calculate_metrics(df_test3)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_pr3, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_pr3, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_pr3, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_pr3, 2)}%')

RMSE:      28133.27       
MAE:       19643.92       
MAPE:      7.76%
Accuracy:  92.24%

# Create a column of the daily error
df_train3['Error'] = df_train3['yhat'] - df_train3['y']
df_train3['Abs_Error'] = np.abs(df_train3['yhat'] - df_train3['y'])

# Create a temporary dataset
df_temp = df_train3[['ds','y','yhat','Error','Abs_Error']]

# Display the 30 largest prediction errors
df_temp.nlargest(30,'Abs_Error')

# Flag the first day of the month for all months except January if it is not a weekend
dfg1['is_1st_day_non_weekend_excl_jan'] = (
    (dfg1['dayofmonth'] == 1) & 
    (dfg1['dayofweek'] < 5) & 
    (dfg1['month'] != 1)
).astype(int)

# Flag 28 December
dfg1['is_28_december'] = ((dfg1['dayofmonth'] == 28) & (dfg1['month'] == 12)).astype(int)

# Make a copy of df grocery (dfg1) to df predict (dfp)
dfp = dfg1.copy()

# Display the updated DataFrame
dfp.tail()

# Step 1: Drop the existing 'date' column
# The index is the date. The date column should be dropped, to avoid duplication.
if 'date' in dfp.columns:
    dfp = dfp.drop(columns=['date'])

# Step 2: Reset the index to create the 'date' column
dfp = dfp.reset_index()

# Step 3: Rename the columns: 'sales' to 'y' and 'index' (former date index) to 'ds' as is required for Prophet
dfp = dfp.rename(columns={'sales': 'y', 'date': 'ds'})

# Step 4: Keep only the 'ds', 'y' and the chosen regressor columns
dfp = dfp[['ds', 'y', 'first_weekend_of_month', 'first_weekend_after_15th', 'dcoilwtico', 
           'onpromotion', 'log_onpromotion', 'earthquake_1_day_after', 'eq_2_days_after', 
            'eq_3_days_after', 'eq_4_days_after', 'eq_5_days_after', 'eq_6_days_after', 
           'is_1st_day_non_weekend_excl_jan', 'is_28_december']]

# Step 5: Ensure 'y' is numeric
dfp['y'] = pd.to_numeric(dfp['y'], errors='coerce')

df_train32 = dfp[(dfp['ds'] < test_start_date) ]
#df_test = df_pred[df_pred['ds'] >= test_start_date]

df_train32.head()

# Create an instance of the model and add the holiday data to it
pr32 = Prophet(holidays=ecuador_holidays_df,
              seasonality_mode='additive')

# Add the regressors to the model
pr32.add_regressor('first_weekend_of_month', mode='multiplicative')
pr32.add_regressor('first_weekend_after_15th', mode='multiplicative')
pr32.add_regressor('dcoilwtico')
#pr32.add_regressor('onpromotion')
#pr32.add_regressor('log_onpromotion')
pr32.add_regressor('earthquake_1_day_after')
pr32.add_regressor('eq_2_days_after')
pr32.add_regressor('eq_3_days_after')
pr32.add_regressor('eq_4_days_after')
pr32.add_regressor('eq_5_days_after')
pr32.add_regressor('eq_6_days_after')
pr32.add_regressor('is_1st_day_non_weekend_excl_jan')
pr32.add_regressor('is_28_december')

<prophet.forecaster.Prophet at 0x1c0904e4550>

# Train the model
pr32.fit(df_train32)

13:11:26 - cmdstanpy - INFO - Chain [1] start processing
13:11:26 - cmdstanpy - INFO - Chain [1] done processing

<prophet.forecaster.Prophet at 0x1c0904e4550>

# Make a future data frame for the number of specified days. 
# This is seven and a half months
future32 = pr32.make_future_dataframe(periods=num_days)

# Merge the regressor columns onto the future dataframe
# Ensure that the values for 'first_weekend_of_month' and 'first_weekend_after_15th' 
# are correctly populated for the future dates
# Note: if they are commented out above, they will not be active even though they are included here.
future32 = future32.merge(dfp[['ds', 'first_weekend_of_month', 'first_weekend_after_15th', 'dcoilwtico', 
                             'onpromotion', 'log_onpromotion', 'earthquake_1_day_after', 'eq_2_days_after', 
                             'eq_3_days_after', 'eq_4_days_after', 'eq_5_days_after', 'eq_6_days_after', 
                            'is_1st_day_non_weekend_excl_jan', 'is_28_december']], 
                        on='ds', 
                        how='left')

# Predict sales for the next two months
df_pred32 = pr32.predict(future32)

plt.figure(figsize=(14, 8))  # Set the plot size to 14x8
fig1 = pr32.plot(df_pred32)
plt.show()

<Figure size 1400x800 with 0 Axes>

plt.figure(figsize=(14, 8))
pr32.plot_components(df_pred32)
plt.show()

<Figure size 1400x800 with 0 Axes>

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# Filter out the test set.
df_test32 = df_pred32[df_pred32['ds'] >= test_start_date]

# Confirm the test set start date
test_start_date.strftime("%Y-%m-%d")

'2017-01-01'

# Ensure df is a standalone copy
df_test32 = df_test32.copy()

# Convert 'ds' to datetime in both dataframes to ensure consistency
df_test32['ds'] = pd.to_datetime(df_test32['ds'])
dfp['ds'] = pd.to_datetime(dfp['ds'])

df_test32 = df_test32.merge(dfp[['ds', 'y']], on='ds', how='left')

# Merge predictions with actual data over the entire dataset
df_pred32_long = df_pred32.merge(dfp[['ds', 'y']], on='ds', how='left')

# Create a plot dataset showing the last two years of the predictions.
# This comprises both the train and test set.
df_pred32_plt = df_pred32_long[df_pred32_long['ds'] >= '2015-08-15']

model_name = 'Pr3.2'
plot_train_set(test_start_date, df_pred32_plt, model_name)

df_train32.tail()

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# Filter for training data only, before the cutoff date of 2017-01-01
df_train32 = df_pred32[df_pred32['ds'] < test_start_date]

# Merge the actual sales values to the predicted test set
df_train32 = df_pred32[['ds', 'yhat']].merge(dfp, on='ds', how='left')

# Call the function to calculate the performance metrics
rmse_trn32, mae_trn32, mape_trn32, accuracy_trn32 = calculate_metrics(df_train32)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_trn32, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_trn32, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_trn32, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_trn32, 2)}%')

RMSE:      22232.24       
MAE:       15765.8        
MAPE:      17.08%
Accuracy:  82.92%

df_test32 = model_new_years_day_prophet(df_train, df_test32)

13:11:29 - cmdstanpy - INFO - Chain [1] start processing
13:11:29 - cmdstanpy - INFO - Chain [1] done processing

plot_test_set(model_name, df_pred32_plt, df_test32)

# Call the function to calculate the performance metrics
rmse_pr32, mae_pr32, mape_pr32, accuracy_pr32 = calculate_metrics(df_test32)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_pr32, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_pr32, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_pr32, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_pr32, 2)}%')

RMSE:      26952.28       
MAE:       18853.49       
MAPE:      7.49%
Accuracy:  92.51%

# Make a copy of df grocery (dfg) to df predict (dfp)
dfp = dfg1.copy()

dfp.tail()

# Step 1: Drop the existing 'date' column
# The index is the date. The date column should be dropped, to avoid duplication.
if 'date' in dfp.columns:
    dfp = dfp.drop(columns=['date'])

# Step 2: Reset the index to create the 'date' column
dfp = dfp.reset_index()

# Step 3: Rename the columns: 'sales' to 'y' and 'index' (former date index) to 'ds' as is required for Prophet
dfp = dfp.rename(columns={'sales': 'y', 'date': 'ds'})

# Step 4: Keep only the 'ds', 'y' and the chosen regressor columns
dfp = dfp[['ds', 'y', 'first_weekend_of_month', 'first_weekend_after_15th', 'dcoilwtico', 
           'onpromotion', 'log_onpromotion', 'earthquake_1_day_after', 'eq_2_days_after', 
            'eq_3_days_after', 'eq_4_days_after', 'eq_5_days_after', 'eq_6_days_after', 
           'is_1st_day_non_weekend_excl_jan', 'is_28_december']]

# Ensure 'y' is numeric
dfp['y'] = pd.to_numeric(dfp['y'], errors='coerce')

# Display the modified dataframe's info
dfp.info()

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 1684 entries, 0 to 1683
Data columns (total 15 columns):
 #   Column                           Non-Null Count  Dtype         
---  ------                           --------------  -----         
 0   ds                               1684 non-null   datetime64[ns]
 1   y                                1684 non-null   float64       
 2   first_weekend_of_month           1684 non-null   int64         
 3   first_weekend_after_15th         1684 non-null   int64         
 4   dcoilwtico                       1684 non-null   float64       
 5   onpromotion                      1684 non-null   float64       
 6   log_onpromotion                  1684 non-null   float64       
 7   earthquake_1_day_after           1684 non-null   int64         
 8   eq_2_days_after                  1684 non-null   int64         
 9   eq_3_days_after                  1684 non-null   int64         
 10  eq_4_days_after                  1684 non-null   int64         
 11  eq_5_days_after                  1684 non-null   int64         
 12  eq_6_days_after                  1684 non-null   int64         
 13  is_1st_day_non_weekend_excl_jan  1684 non-null   int32         
 14  is_28_december                   1684 non-null   int32         
dtypes: datetime64[ns](1), float64(4), int32(2), int64(8)
memory usage: 184.3 KB

df_train4 = dfp[(dfp['ds'] < test_start_date) ]
df_test4 = dfp[(dfp['ds'] >= test_start_date) ]

# Make copies to avoid warnings when scaling these data frames.
df_train4 = df_train4.copy()
df_test4 = df_test4.copy()

# Define the list of continuous regressors
continuous_regressors = ['dcoilwtico', 'onpromotion', 'log_onpromotion']

# StandardScaler standardizes features by removing the mean and scaling to unit variance.
scaler = StandardScaler()

# View the header
df_train4.head()

# Apply the scaler to the continuous regressors in the training dataset
# The `fit_transform` method calculates the mean and standard deviation for scaling 
# and applies the transformation to standardize the data.
df_train4[continuous_regressors] = scaler.fit_transform(
    df_train4[continuous_regressors]) # Subset of the DataFrame containing only continuous regressors

# View the heder after scaling
df_train4.head()

# Apply the previously fitted StandardScaler to the test dataset
# The `transform` method uses the mean and standard deviation calculated from the training data 
# (stored in the scaler) to scale the continuous regressors in the test dataset. 
# This ensures that the test data is standardized using the same parameters as the training data,
# maintaining consistency for model evaluation.
df_test4[continuous_regressors] = scaler.transform(df_test4[continuous_regressors])

df_test4.head()

# Concatenate df_train4 and df_test4 into dfp_scaled
dfp_scaled = pd.concat([df_train4, df_test4], ignore_index=True)

# Create an instance of the model and add the holiday data to it
pr4 = Prophet(holidays=ecuador_holidays_df,
              seasonality_mode='additive')

# Add regressors to the model
pr4.add_regressor('first_weekend_of_month', mode='multiplicative')
pr4.add_regressor('first_weekend_after_15th', mode='multiplicative')
pr4.add_regressor('dcoilwtico')  # Scaled
#pr4.add_regressor('onpromotion')  # Scaled
#pr4.add_regressor('log_onpromotion')  # Scaled
pr4.add_regressor('earthquake_1_day_after')
pr4.add_regressor('eq_2_days_after')
pr4.add_regressor('eq_3_days_after')
pr4.add_regressor('eq_4_days_after')
pr4.add_regressor('eq_5_days_after')
pr4.add_regressor('eq_6_days_after')
pr4.add_regressor('is_1st_day_non_weekend_excl_jan')
pr4.add_regressor('is_28_december')

<prophet.forecaster.Prophet at 0x1c097d780a0>

# Train the model
pr4.fit(df_train4)

13:11:30 - cmdstanpy - INFO - Chain [1] start processing
13:11:30 - cmdstanpy - INFO - Chain [1] done processing

<prophet.forecaster.Prophet at 0x1c097d780a0>

# Make a future data frame for the number of specified days. 
# This is seven and a half months
future4 = pr4.make_future_dataframe(periods=num_days)

# Merge the regressor columns onto the future dataframe
# Note: if they are commented out above, they will not be active even though they are included here.
future4 = future4.merge(dfp_scaled[['ds', 'first_weekend_of_month', 'first_weekend_after_15th', 'dcoilwtico',
                             'onpromotion', 'log_onpromotion', 'earthquake_1_day_after', 'eq_2_days_after',
                             'eq_3_days_after', 'eq_4_days_after', 'eq_5_days_after', 'eq_6_days_after',
                             'is_1st_day_non_weekend_excl_jan', 'is_28_december']],
                        on='ds', how='left')

# Predict sales for the next two months
df_pred4 = pr4.predict(future4)

plt.figure(figsize=(14, 8))
fig1 = pr4.plot(df_pred4)
plt.show()

<Figure size 1400x800 with 0 Axes>

plt.figure(figsize=(14, 8))
pr4.plot_components(df_pred4)
plt.show()

<Figure size 1400x800 with 0 Axes>

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# Filter out the test set.
df_test4 = df_pred4[df_pred4['ds'] >= test_start_date]

# Confirm the test set start date
test_start_date.strftime("%Y-%m-%d")

'2017-01-01'

# Ensure df_test4 is a standalone copy
df_test4 = df_test4.copy()

# Convert 'ds' to datetime in both dataframes to ensure consistency
df_test4['ds'] = pd.to_datetime(df_test4['ds'])
dfp['ds'] = pd.to_datetime(dfp['ds'])

df_test4 = df_test4.merge(dfp[['ds', 'y']], on='ds', how='left')

# Merge predictions with actual data over the entire dataset
df_pred4_long = df_pred4.merge(dfp[['ds', 'y']], on='ds', how='left')

# Create a plot dataset showing the last two years of the predictions.
# This comprises both the train and test set.
df_pred4_plt = df_pred4_long[df_pred4_long['ds'] >= '2015-08-15']

model_name = 'Pr4'
plot_train_set(test_start_date, df_pred4_plt, model_name)

df_train4.tail()

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# Filter for training data only, before the cutoff date of 2017-01-01
df_train4 = df_pred4[df_pred4['ds'] < test_start_date]

# Merge the actual sales values to the predicted test set
df_train4 = df_pred4[['ds', 'yhat']].merge(dfp, on='ds', how='left')

# Call the function to calculate the performance metrics
rmse_trn4, mae_trn4, mape_trn4, accuracy_trn4 = calculate_metrics(df_train4)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_trn4, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_trn4, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_trn4, 1)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_trn4, 1)}%')

RMSE:      22232.9        
MAE:       15763.88       
MAPE:      17.1%
Accuracy:  82.9%

df_test4 = model_new_years_day_prophet(df_train, df_test4)

13:11:33 - cmdstanpy - INFO - Chain [1] start processing
13:11:33 - cmdstanpy - INFO - Chain [1] done processing

plot_test_set(model_name, df_pred4_plt, df_test4)

# Call the function to calculate the performance metrics
rmse_pr4, mae_pr4, mape_pr4, accuracy_pr4 = calculate_metrics(df_test4)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_pr4, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_pr4, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_pr4, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_pr4, 2)}%')

RMSE:      26955.31       
MAE:       18852.56       
MAPE:      7.49%
Accuracy:  92.51%

# Make a copy of df grocery (dfg1) to df predict (dfp)
dfp = dfg1.copy()

# Step 1: Drop the existing 'date' column
# The index is the date. The date column should be dropped, to avoid duplication.
if 'date' in dfp.columns:
    dfp = dfp.drop(columns=['date'])

# Step 2: Reset the index to create the 'date' column
dfp = dfp.reset_index()

# Step 3: Rename the columns: 'sales' to 'y' and 'index' (former date index) to 'ds' as is required for Prophet
dfp = dfp.rename(columns={'sales': 'y', 'date': 'ds'})

# Step 4: Keep only the 'ds', 'y' and the chosen regressor columns
dfp = dfp[['ds', 'y', 'first_weekend_of_month', 'first_weekend_after_15th', 'dcoilwtico', 
           'onpromotion', 'log_onpromotion', 'earthquake_1_day_after', 'eq_2_days_after', 
            'eq_3_days_after', 'eq_4_days_after', 'eq_5_days_after', 'eq_6_days_after', 
           'is_1st_day_non_weekend_excl_jan', 'is_28_december']]

# Step 5: Ensure 'y' is numeric
dfp['y'] = pd.to_numeric(dfp['y'], errors='coerce')

df_train5 = dfp[(dfp['ds'] < test_start_date) ]

# Set up an initial parameter grid with an extensive list of parameters to be tuned

param_grid = {
    'changepoint_prior_scale': [0.005, 0.01, 0.05, 0.1, 0.2],  # flexibility of trend changepoints
    'seasonality_prior_scale': [0.1, 1.0, 5.0, 10.0],  # flexibility of seasonal components
    'holidays_prior_scale': [0.05, 0.1, 1.0, 5.0, 10.0],  # flexibility of holidays
    'seasonality_mode': ['additive', 'multiplicative'],  # whether seasonality is additive or multiplicative
    'n_changepoints': [10, 25, 50, 75],  # number of potential trend changepoints
    'yearly_seasonality': [10, 20, 30],  # Fourier order for yearly seasonality
    'weekly_seasonality': [2, 3, 5, 7]  # Fourier order for weekly seasonality
}


# Generate all combinations of parameters
all_params = [dict(zip(param_grid.keys(), v)) for v in itertools.product(*param_grid.values())]
rmses = []  # Store the RMSEs for all parameter cominations here

"""
# Record the start time
start_time = time.time()

# Use cross validation to evaluate all parameters
for params in all_params:
    # Create an instance of the model and add the holiday data to it
    pr5 = Prophet(**params, holidays=ecuador_holidays_df)
    
    # Add the regressors to the model
    pr5.add_regressor('first_weekend_of_month', mode='multiplicative')
    pr5.add_regressor('first_weekend_after_15th', mode='multiplicative')
    pr5.add_regressor('dcoilwtico')
    #pr5.add_regressor('onpromotion'')
    #pr5.add_regressor('log_onpromotion')
    pr5.add_regressor('earthquake_1_day_after')
    pr5.add_regressor('eq_2_days_after')
    pr5.add_regressor('eq_3_days_after')
    pr5.add_regressor('eq_4_days_after')
    pr5.add_regressor('eq_5_days_after')
    pr5.add_regressor('eq_6_days_after')
    pr5.add_regressor('is_1st_day_non_weekend_excl_jan')
    pr5.add_regressor('is_28_december')
    
    pr5.fit(df_train5)  # Fit model with given params
    
    df_cv = cross_validation(pr5, initial='1095 days', period='31 days', horizon = '227 days', parallel="processes")
    df_p = performance_metrics(df_cv, rolling_window=1)
    rmses.append(df_p['rmse'].values[0])

# Record the end time
end_time = time.time()

# Calculate and print the total execution time
execution_time = end_time - start_time

# Convert to hours, minutes, and seconds
hours, remainder = divmod(execution_time, 3600)
minutes, seconds = divmod(remainder, 60)

# Print the formatted execution time
print(f"Total execution time: {int(hours)} hours, {int(minutes)} minutes, and {seconds:.2f} seconds")
"""

'\n# Record the start time\nstart_time = time.time()\n\n# Use cross validation to evaluate all parameters\nfor params in all_params:\n    # Create an instance of the model and add the holiday data to it\n    pr5 = Prophet(**params, holidays=ecuador_holidays_df)\n    \n    # Add the regressors to the model\n    pr5.add_regressor(\'first_weekend_of_month\', mode=\'multiplicative\')\n    pr5.add_regressor(\'first_weekend_after_15th\', mode=\'multiplicative\')\n    pr5.add_regressor(\'dcoilwtico\')\n    #pr5.add_regressor(\'onpromotion\'\')\n    #pr5.add_regressor(\'log_onpromotion\')\n    pr5.add_regressor(\'earthquake_1_day_after\')\n    pr5.add_regressor(\'eq_2_days_after\')\n    pr5.add_regressor(\'eq_3_days_after\')\n    pr5.add_regressor(\'eq_4_days_after\')\n    pr5.add_regressor(\'eq_5_days_after\')\n    pr5.add_regressor(\'eq_6_days_after\')\n    pr5.add_regressor(\'is_1st_day_non_weekend_excl_jan\')\n    pr5.add_regressor(\'is_28_december\')\n    \n    pr5.fit(df_train5)  # Fit model with given params\n    \n    df_cv = cross_validation(pr5, initial=\'1095 days\', period=\'31 days\', horizon = \'227 days\', parallel="processes")\n    df_p = performance_metrics(df_cv, rolling_window=1)\n    rmses.append(df_p[\'rmse\'].values[0])\n\n# Record the end time\nend_time = time.time()\n\n# Calculate and print the total execution time\nexecution_time = end_time - start_time\n\n# Convert to hours, minutes, and seconds\nhours, remainder = divmod(execution_time, 3600)\nminutes, seconds = divmod(remainder, 60)\n\n# Print the formatted execution time\nprint(f"Total execution time: {int(hours)} hours, {int(minutes)} minutes, and {seconds:.2f} seconds")\n'

"""
# Find the best parameters
tuning_results = pd.DataFrame(all_params)
tuning_results['rmse'] = rmses
"""

"\n# Find the best parameters\ntuning_results = pd.DataFrame(all_params)\ntuning_results['rmse'] = rmses\n"

"""
# Save the DataFrame as a CSV file
tuning_results.to_csv('./data/prophet_5_grid_search_results.csv', index=False)

# Save the cross-validations results as a CSV file.
df_cv.to_csv('./data/cross_validation_results_pr5.csv', index=False)
"""

"\n# Save the DataFrame as a CSV file\ntuning_results.to_csv('./data/prophet_5_grid_search_results.csv', index=False)\n\n# Save the cross-validations results as a CSV file.\ndf_cv.to_csv('./data/cross_validation_results_pr5.csv', index=False)\n"

# Load the saved CSV files into their respective DataFrames
tuning_results = pd.read_csv('./data/prophet_5_grid_search_results.csv')
df_cv = pd.read_csv('./data/cross_validation_results_pr5.csv')

# Display the top 20 results
tuning_results.sort_values('rmse').head(20)

# Display the best hyperparameters
tuning_results.sort_values('rmse').reset_index(drop=True).iloc[0]

changepoint_prior_scale              0.01
seasonality_prior_scale               0.1
holidays_prior_scale                 0.05
seasonality_mode           multiplicative
n_changepoints                         10
yearly_seasonality                     30
weekly_seasonality                      7
rmse                         41878.212478
Name: 0, dtype: object

# Store the best hyperparameters in a dictionary
best_parameters = dict(tuning_results.sort_values('rmse').reset_index(drop=True).drop('rmse',axis='columns').iloc[0])

best_parameters

{'changepoint_prior_scale': 0.01,
 'seasonality_prior_scale': 0.1,
 'holidays_prior_scale': 0.05,
 'seasonality_mode': 'multiplicative',
 'n_changepoints': 10,
 'yearly_seasonality': 30,
 'weekly_seasonality': 7}

# Convert the saved Columns to Datetime
df_cv['ds'] = pd.to_datetime(df_cv['ds'])
df_cv['cutoff'] = pd.to_datetime(df_cv['cutoff'])

df_p = performance_metrics(df_cv)
df_p

fig = plot_cross_validation_metric(df_cv, metric='rmse')

# Create an instance of the model with the best hyperparameters and add the holiday data to it
pr5 = Prophet(
    holidays=ecuador_holidays_df,
    changepoint_prior_scale=best_parameters['changepoint_prior_scale'], 
    seasonality_prior_scale=best_parameters['seasonality_prior_scale'],
    holidays_prior_scale=best_parameters.get('holidays_prior_scale'),
    seasonality_mode=best_parameters.get('seasonality_mode'),
    n_changepoints=best_parameters.get('n_changepoints'),
    yearly_seasonality=best_parameters.get('yearly_seasonality'),  
    weekly_seasonality=best_parameters.get('weekly_seasonality')
)

# Add the regressors to the model
pr5.add_regressor('first_weekend_of_month', mode='multiplicative')
pr5.add_regressor('first_weekend_after_15th', mode='multiplicative')
pr5.add_regressor('dcoilwtico')
#pr5.add_regressor('onpromotion')
#pr5.add_regressor('log_onpromotion')
pr5.add_regressor('earthquake_1_day_after')
pr5.add_regressor('eq_2_days_after')
pr5.add_regressor('eq_3_days_after')
pr5.add_regressor('eq_4_days_after')
pr5.add_regressor('eq_5_days_after')
pr5.add_regressor('eq_6_days_after')
pr5.add_regressor('is_1st_day_non_weekend_excl_jan')
pr5.add_regressor('is_28_december')

<prophet.forecaster.Prophet at 0x1c097516a10>

# Train the model
pr5.fit(df_train5)

13:11:35 - cmdstanpy - INFO - Chain [1] start processing
13:11:36 - cmdstanpy - INFO - Chain [1] done processing

<prophet.forecaster.Prophet at 0x1c097516a10>

# Make a future data frame for the number of specified days. 
# This is seven and a half months
future5 = pr5.make_future_dataframe(periods=num_days)

# Merge the regressor columns onto the future dataframe
# Note: if they are commented out above, they will not be active even though they are included here.
future5 = future5.merge(dfp[['ds', 'first_weekend_of_month', 'first_weekend_after_15th', 'dcoilwtico', 
                             'onpromotion', 'log_onpromotion', 'earthquake_1_day_after', 'eq_2_days_after', 
                             'eq_3_days_after', 'eq_4_days_after', 'eq_5_days_after', 'eq_6_days_after', 
                            'is_1st_day_non_weekend_excl_jan', 'is_28_december']], 
                        on='ds', 
                        how='left')

# Predict sales for the next two months
df_pred5 = pr5.predict(future5)

plt.figure(figsize=(14, 8)) 
fig1 = pr5.plot(df_pred5)
plt.show()

<Figure size 1400x800 with 0 Axes>

plt.figure(figsize=(14, 8))
pr5.plot_components(df_pred5)
plt.show()

<Figure size 1400x800 with 0 Axes>

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# Filter out the test set.
df_test5 = df_pred5[df_pred5['ds'] >= test_start_date]

# Convert 'ds' to datetime in both dataframes to ensure consistency
df_test5.loc[:, 'ds'] = pd.to_datetime(df_test5['ds'])
dfp['ds'] = pd.to_datetime(dfp['ds'])

df_test5 = df_test5.merge(dfp[['ds', 'y']], on='ds', how='left')

# Merge predictions with actual data over the entire dataset
df_pred5_long = df_pred5.merge(dfp[['ds', 'y']], on='ds', how='left')

# Create a plot dataset showing the last two years of the predictions.
# This comprises both the train and test set.
df_pred5_plt = df_pred5_long[df_pred5_long['ds'] >= '2015-08-15']

model_name = 'Pr5'
plot_train_set(test_start_date, df_pred5_plt, model_name)

plot_train_set(test_start_date, df_pred32_plt, 'Pr3.2')

df_train5.tail()

# Prophet has predicted for the entire training set plus the additional 227 days in the test set.
# Filter for training data only, before the cutoff date of 2017-01-01
df_train5 = df_pred5[df_pred5['ds'] < test_start_date]

# Merge the actual sales values to the predicted test set
df_train5 = df_pred5[['ds', 'yhat']].merge(dfp, on='ds', how='left')

# Call the function to calculate the performance metrics
rmse_trn5, mae_trn5, mape_trn5, accuracy_trn5 = calculate_metrics(df_train5)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_trn5, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_trn5, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_trn5, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_trn5, 2)}%')

RMSE:      18190.09       
MAE:       12449.73       
MAPE:      10.61%
Accuracy:  89.39%

df_test5 = model_new_years_day_prophet(df_train, df_test5)

13:11:38 - cmdstanpy - INFO - Chain [1] start processing
13:11:38 - cmdstanpy - INFO - Chain [1] done processing

plot_test_set('Pr5', df_pred5_plt, df_test5)

plot_test_set('Pr3.2', df_pred32_plt, df_test32)

# Call the function to calculate the performance metrics
rmse_pr5, mae_pr5, mape_pr5, accuracy_pr5 = calculate_metrics(df_test5)

# Print the performance metrics
print(f'{"RMSE:":<{label_width}} {round(rmse_pr5, 2):<{value_width}}')
print(f'{"MAE:":<{label_width}} {round(mae_pr5, 2):<{value_width}}')
print(f'{"MAPE:":<{label_width}} {round(mape_pr5, 2)}%')
print(f'{"Accuracy:":<{label_width}} {round(100 - mape_pr5, 2)}%')

RMSE:      24679.15       
MAE:       16990.55       
MAPE:      7.03%
Accuracy:  92.97%

	date	type	locale	locale_name	description	transferred
0	2012-03-02	Holiday	Local	Manta	Fundacion de Manta	False
1	2012-04-01	Holiday	Regional	Cotopaxi	Provincializacion de Cotopaxi	False
2	2012-04-12	Holiday	Local	Cuenca	Fundacion de Cuenca	False
3	2012-04-14	Holiday	Local	Libertad	Cantonizacion de Libertad	False
4	2012-04-21	Holiday	Local	Riobamba	Cantonizacion de Riobamba	False

	date	dcoilwtico
0	1/01/2013	NaN
1	2/01/2013	93.14
2	3/01/2013	92.97
3	4/01/2013	93.12
4	7/01/2013	93.20

	year_month	store_nbr	family	sales	onpromotion
99787	2017-08-01	54	POULTRY	991.986996	17
99788	2017-08-01	54	PREPARED FOODS	1340.000000	0
99789	2017-08-01	54	PRODUCE	10351.598000	263
99790	2017-08-01	54	SCHOOL AND OFFICE SUPPLIES	0.000000	0
99791	2017-08-01	54	SEAFOOD	38.000000	2

family	AUTOMOTIVE	BABY CARE	BEAUTY	BEVERAGES	BOOKS	BREAD/BAKERY	CELEBRATION	CLEANING	DAIRY	DELI	EGGS	FROZEN FOODS	GROCERY I	GROCERY II	HARDWARE	HOME AND KITCHEN I	HOME AND KITCHEN II	HOME APPLIANCES	HOME CARE	LADIESWEAR	LAWN AND GARDEN	LINGERIE	LIQUOR,WINE,BEER	MAGAZINES	MEATS	PERSONAL CARE	PET SUPPLIES	PLAYERS AND ELECTRONICS	POULTRY	PREPARED FOODS	PRODUCE	SCHOOL AND OFFICE SUPPLIES	SEAFOOD
date
2013-01-01	0.0	0.0	2.0	810.0	0.0	180.58900	0.0	186.0	143.0	71.090000	46.0	29.654999	700.000	15.0	0.0	0.0	0.0	0.0	0.0	0.0	2.0	5.0	105.0	0.0	110.801000	25.0	0.0	0.0	42.637000	37.847000	0.00000	0.0	0.000000
2013-01-02	255.0	0.0	207.0	72092.0	0.0	26246.31900	0.0	74629.0	23381.0	15754.500000	10932.0	7115.388000	202020.000	1476.0	46.0	0.0	0.0	21.0	0.0	0.0	121.0	469.0	2411.0	0.0	20871.464028	17204.0	0.0	0.0	13975.884938	5338.111976	0.00000	0.0	1526.750002
2013-01-03	161.0	0.0	125.0	52105.0	0.0	18456.48002	0.0	55893.0	18001.0	11172.455000	7358.0	4760.805009	144878.000	1048.0	37.0	0.0	0.0	15.0	0.0	0.0	83.0	366.0	2476.0	0.0	16597.398113	12568.0	0.0	0.0	10674.393983	3591.388005	0.00000	0.0	1094.310994
2013-01-04	169.0	0.0	133.0	54167.0	0.0	16721.96901	0.0	52064.0	18148.0	10143.209000	6760.0	4525.934004	135754.000	1031.0	57.0	0.0	0.0	13.0	0.0	0.0	127.0	382.0	4796.0	0.0	21625.963055	11303.0	0.0	0.0	10772.515038	4472.965990	0.00000	0.0	1293.120995
2013-01-05	342.0	0.0	191.0	77818.0	0.0	22367.76108	0.0	70128.0	23082.0	13734.945010	8576.0	5781.605010	188356.000	1273.0	87.0	0.0	0.0	11.0	0.0	0.0	180.0	458.0	6715.0	0.0	20879.091050	16819.0	0.0	0.0	13475.009055	5830.073020	0.00000	0.0	1245.637004
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
2017-08-11	441.0	7.0	351.0	189111.0	0.0	28204.48600	870.0	55682.0	44909.0	18680.488004	9428.0	9473.183000	238221.922	1653.0	83.0	2189.0	1504.0	10.0	13747.0	460.0	1156.0	445.0	11033.0	500.0	26944.092010	15978.0	587.0	654.0	25318.297990	5199.494021	118738.14300	3523.0	1272.615997
2017-08-12	403.0	3.0	369.0	182318.0	0.0	27403.54800	849.0	57935.0	42854.0	15689.893010	9560.0	9077.454000	238627.633	1555.0	89.0	2046.0	1504.0	8.0	14692.0	599.0	978.0	462.0	10058.0	483.0	17775.650112	15903.0	508.0	712.0	19134.510058	4573.465992	111788.35090	3644.0	1028.030006
2017-08-13	481.0	9.0	433.0	202354.0	1.0	33262.75700	544.0	61757.0	50112.0	17593.274020	11478.0	6476.555990	258442.452	1504.0	108.0	1874.0	1625.0	2.0	16564.0	722.0	750.0	470.0	4294.0	469.0	18633.015039	18188.0	541.0	741.0	20509.265004	4941.509018	125923.80240	3718.0	1118.047000
2017-08-14	292.0	9.0	337.0	174832.0	0.0	28002.45700	594.0	55032.0	43050.0	14662.741000	14600.0	5422.669000	225355.694	1371.0	66.0	1682.0	1668.0	4.0	13818.0	487.0	654.0	418.0	2839.0	457.0	16975.686040	15077.0	417.0	500.0	18597.508060	4647.375002	115257.59598	2826.0	970.679999
2017-08-15	337.0	8.0	339.0	170773.0	0.0	29158.19502	597.0	58474.0	40707.0	14935.453000	9624.0	5309.611004	224208.125	1300.0	57.0	1669.0	1448.0	3.0	13644.0	503.0	679.0	458.0	3461.0	461.0	17928.170944	14787.0	364.0	592.0	17586.709986	4641.522980	125108.97100	2530.0	970.177005

	sales	onpromotion	dcoilwtico	log_onpromotion
count	1684.000000	1684.000000	1684.000000	1684.000000
mean	203956.493400	21.056578	67.894634	2.044123
std	61832.884741	27.770375	25.690236	1.606593
min	700.000000	0.000000	26.190000	0.000000
25%	157759.750000	0.000000	46.378333	0.000000
50%	194681.310000	6.361111	53.285000	1.996210
75%	239845.256250	39.185185	95.785000	3.693497
max	583222.480000	155.777778	110.620000	5.054829

	store_nbr	city	state	type	cluster
0	1	Quito	Pichincha	D	13
1	2	Quito	Pichincha	D	13
2	3	Quito	Pichincha	D	8
3	4	Quito	Pichincha	D	9
4	5	Santo Domingo	Santo Domingo de los Tsachilas	D	4

	Total	Percent	Types
id	0	0.0	int64
date	0	0.0	object
store_nbr	0	0.0	int64
family	0	0.0	object
sales	0	0.0	float64
onpromotion	0	0.0	int64

	Total	Percent	Types
date	0	0.000000	datetime64[ns]
sales	0	0.000000	float64
onpromotion	0	0.000000	float64
dcoilwtico	481	28.562945	float64

	sales	onpromotion	dcoilwtico	log_onpromotion	date	hour	dayofweek	weekday	quarter	month	month_name	year	dayofyear	dayofmonth	weekofyear	date_offset	season	first_weekend_of_month	first_weekend_after_15th
date
2013-01-01	700.0	0.0	93.140000	0.0	2013-01-01	0	1	Tuesday	1	1	January	2013	1	1	1	1081	Wet	0	0
2013-01-02	202020.0	0.0	93.140000	0.0	2013-01-02	0	2	Wednesday	1	1	January	2013	2	2	1	1082	Wet	0	0
2013-01-03	144878.0	0.0	92.970000	0.0	2013-01-03	0	3	Thursday	1	1	January	2013	3	3	1	1083	Wet	0	0
2013-01-04	135754.0	0.0	93.120000	0.0	2013-01-04	0	4	Friday	1	1	January	2013	4	4	1	1084	Wet	0	0
2013-01-05	188356.0	0.0	93.146667	0.0	2013-01-05	0	5	Saturday	1	1	January	2013	5	5	1	1085	Wet	1	0
2013-01-06	207297.0	0.0	93.173333	0.0	2013-01-06	0	6	Sunday	1	1	January	2013	6	6	1	1086	Wet	1	0

	ds	y
1452	2016-12-27	242792.095
1453	2016-12-28	265759.513
1454	2016-12-29	257006.744
1455	2016-12-30	329599.654
1456	2016-12-31	302232.811

	ds	trend	yhat_lower	yhat_upper	trend_lower	trend_upper	additive_terms	additive_terms_lower	additive_terms_upper	weekly	weekly_lower	weekly_upper	yearly	yearly_lower	yearly_upper	yhat
0	2013-01-01	147120.786692	115324.889066	194569.477068	147120.786692	147120.786692	6465.581964	6465.581964	6465.581964	-19223.831061	-19223.831061	-19223.831061	25689.413025	25689.413025	25689.413025	153586.368656
1	2013-01-02	147196.805349	107531.416442	187407.456971	147196.805349	147196.805349	-63.132652	-63.132652	-63.132652	-19759.400670	-19759.400670	-19759.400670	19696.268018	19696.268018	19696.268018	147133.672697
2	2013-01-03	147272.824005	82178.953317	164483.587504	147272.824005	147272.824005	-24320.567781	-24320.567781	-24320.567781	-38164.714351	-38164.714351	-38164.714351	13844.146569	13844.146569	13844.146569	122952.256223
3	2013-01-04	147348.842661	91646.242643	172246.654461	147348.842661	147348.842661	-16196.996936	-16196.996936	-16196.996936	-24411.660408	-24411.660408	-24411.660408	8214.663472	8214.663472	8214.663472	131151.845725
4	2013-01-05	147424.861317	153286.333676	231219.211937	147424.861317	147424.861317	43784.377048	43784.377048	43784.377048	40901.410256	40901.410256	40901.410256	2882.966793	2882.966793	2882.966793	191209.238365

	ds	trend	yhat_lower	yhat_upper	trend_lower	trend_upper	additive_terms	additive_terms_lower	additive_terms_upper	weekly	weekly_lower	weekly_upper	yearly	yearly_lower	yearly_upper	yhat	y
0	2017-01-01	244974.255743	290760.384698	372051.725858	244974.255743	244974.255743	87858.666352	87858.666352	87858.666352	62169.253327	62169.253327	62169.253327	25689.413025	25689.413025	25689.413025	332832.922095	1998.000
1	2017-01-02	245021.669456	221004.414431	303284.180416	245021.669456	245021.669456	18185.210925	18185.210925	18185.210925	-1511.057093	-1511.057093	-1511.057093	19696.268018	19696.268018	19696.268018	263206.880382	415182.950
2	2017-01-03	245069.083170	200397.990576	279405.872380	245069.083170	245069.083170	-5379.684492	-5379.684492	-5379.684492	-19223.831061	-19223.831061	-19223.831061	13844.146569	13844.146569	13844.146569	239689.398678	321085.262
3	2017-01-04	245116.496883	194725.344188	272302.742938	245116.496883	245116.496883	-11544.737198	-11544.737198	-11544.737198	-19759.400670	-19759.400670	-19759.400670	8214.663472	8214.663472	8214.663472	233571.759685	274170.595
4	2017-01-05	245163.910596	166589.014252	251780.696684	245163.910596	245163.910596	-35281.747558	-35281.747558	-35281.747558	-38164.714351	-38164.714351	-38164.714351	2882.966793	2882.966793	2882.966793	209882.163038	226863.181

	holiday	ds	lower_window	upper_window
25	Año Nuevo	2013-01-01	0	0
26	Carnaval	2013-02-11	-1	2
27	Carnaval	2013-02-12	-1	2
28	Viernes Santo	2013-03-29	-1	2
29	Día del Trabajo	2013-05-01	-1	2
30	Batalla de Pichincha	2013-05-24	-1	2
31	Primer Grito de Independencia	2013-08-10	-1	2
32	Independencia de Guayaquil	2013-10-09	-1	2
33	Día de los Difuntos	2013-11-02	-1	2
34	Independencia de Cuenca	2013-11-03	-1	2
35	Día de Navidad	2013-12-25	-1	2
36	Año Nuevo	2014-01-01	0	0
37	Carnaval	2014-03-03	-1	2
38	Carnaval	2014-03-04	-1	2
39	Viernes Santo	2014-04-18	-1	2
40	Día del Trabajo	2014-05-01	-1	2
41	Batalla de Pichincha	2014-05-24	-1	2
42	Primer Grito de Independencia	2014-08-10	-1	2
43	Independencia de Guayaquil	2014-10-09	-1	2
44	Día de los Difuntos	2014-11-02	-1	2
45	Independencia de Cuenca	2014-11-03	-1	2
46	Día de Navidad	2014-12-25	-1	2
47	Año Nuevo	2015-01-01	0	0
48	Carnaval	2015-02-16	-1	2
49	Carnaval	2015-02-17	-1	2

	ds	trend	yhat_lower	yhat_upper	trend_lower	trend_upper	Año Nuevo	Año Nuevo_lower	Año Nuevo_upper	additive_terms	additive_terms_lower	additive_terms_upper	holidays	holidays_lower	holidays_upper	weekly	weekly_lower	weekly_upper	yearly	yearly_lower	yearly_upper	yhat
0	2013-01-01	146149.138041	-66922.134382	3201.276334	146149.138041	146149.138041	-201583.851542	-201583.851542	-201583.851542	-177985.424157	-177985.424157	-177985.424157	-201583.851542	-201583.851542	-201583.851542	-19373.952292	-19373.952292	-19373.952292	42972.379677	42972.379677	42972.379677	-31836.286116
1	2013-01-02	146238.742185	125011.658248	194727.875682	146238.742185	146238.742185	0.000000	0.000000	0.000000	16123.223823	16123.223823	16123.223823	0.000000	0.000000	0.000000	-20348.514956	-20348.514956	-20348.514956	36471.738779	36471.738779	36471.738779	162361.966008
2	2013-01-03	146328.346330	100065.329218	170258.540955	146328.346330	146328.346330	0.000000	0.000000	0.000000	-9768.256527	-9768.256527	-9768.256527	0.000000	0.000000	0.000000	-39674.946275	-39674.946275	-39674.946275	29906.689748	29906.689748	29906.689748	136560.089803
3	2013-01-04	146417.950475	110660.366055	181659.410299	146417.950475	146417.950475	0.000000	0.000000	0.000000	-650.600700	-650.600700	-650.600700	0.000000	0.000000	0.000000	-24027.011310	-24027.011310	-24027.011310	23376.410610	23376.410610	23376.410610	145767.349774
4	2013-01-05	146507.554620	167538.561910	240038.433615	146507.554620	146507.554620	0.000000	0.000000	0.000000	57917.918814	57917.918814	57917.918814	0.000000	0.000000	0.000000	40940.274164	40940.274164	40940.274164	16977.644649	16977.644649	16977.644649	204425.473433

	ds	y	yhat	Error	Abs_Error
90	2017-04-01	474003.906	283279.328716	-190724.577284	190724.577284
0	2017-01-01	1998.000	149224.332047	147226.332047	147226.332047
1	2017-01-02	415182.950	279197.893568	-135985.056432	135985.056432
154	2017-06-04	420484.196	314419.044956	-106065.151044	106065.151044
182	2017-07-02	398356.868	308805.184488	-89551.683512	89551.683512

	sales	onpromotion	dcoilwtico	log_onpromotion
date
2017-08-11	238221.922	43.500000	48.810000	3.795489
2017-08-12	238627.633	43.055556	48.403333	3.785451
2017-08-13	258442.452	42.962963	47.996667	3.783348
2017-08-14	225355.694	43.240741	47.590000	3.789646
2017-08-15	224208.125	42.518519	47.570000	3.773187

	ds	y	yhat	Error	Abs_Error
0	2017-01-01	1998.000	332832.922095	330834.922095	330834.922095
90	2017-04-01	474003.906	284095.517660	-189908.388340	189908.388340
1	2017-01-02	415182.950	263206.880382	-151976.069618	151976.069618
120	2017-05-01	400850.633	253164.213753	-147686.419247	147686.419247
154	2017-06-04	420484.196	311331.178890	-109153.017110	109153.017110

	ds	y	yhat	Error	Abs_Error
90	2017-04-01	474003.906	283279.328716	-190724.577284	190724.577284
1	2017-01-02	415182.950	279197.893568	-135985.056432	135985.056432
154	2017-06-04	420484.196	314419.044956	-106065.151044	106065.151044
182	2017-07-02	398356.868	308805.184488	-89551.683512	89551.683512
181	2017-07-01	375440.305	286122.290548	-89318.014452	89318.014452
145	2017-05-26	284289.827	196917.953185	-87371.873815	87371.873815
126	2017-05-07	376379.473	298316.475723	-78062.997277	78062.997277
212	2017-08-01	301453.492	224688.817338	-76764.674662	76764.674662
153	2017-06-03	366851.135	292745.704174	-74105.430826	74105.430826
91	2017-04-02	380816.880	308252.939783	-72563.940217	72563.940217

	holiday	ds	lower_window	upper_window
1	Carnaval	2016-02-08	-1	2
2	Carnaval	2016-02-09	-1	2
3	Viernes Santo	2016-03-25	-1	2
4	Día del Trabajo	2016-05-01	-1	2
5	Batalla de Pichincha	2016-05-24	-1	2
6	Primer Grito de Independencia	2016-08-10	-1	2
7	Independencia de Guayaquil	2016-10-09	-1	2
8	Día de los Difuntos	2016-11-02	-1	2
9	Independencia de Cuenca	2016-11-03	-1	2
12	Año Nuevo (observado)	2017-01-02	-1	2
13	Carnaval	2017-02-27	-1	2
14	Carnaval	2017-02-28	-1	2
15	Viernes Santo	2017-04-14	-1	2
16	Día del Trabajo	2017-05-01	-1	2
17	Batalla de Pichincha	2017-05-24	-1	2

	ds	trend	yhat_lower	yhat_upper	trend_lower	trend_upper	Año Nuevo	Año Nuevo_lower	Año Nuevo_upper	additive_terms	additive_terms_lower	additive_terms_upper	extra_regressors_multiplicative	extra_regressors_multiplicative_lower	extra_regressors_multiplicative_upper	multiplicative_terms	multiplicative_terms_lower	multiplicative_terms_upper	weekly	weekly_lower	weekly_upper	yearly	yearly_lower	yearly_upper	yhat
0	2013-01-01	147105.039685	-26027.943476	43924.832612	147105.039685	147105.039685	-1.08038	-1.08038	-1.08038	19521.127696	19521.127696	19521.127696	-1.08038	-1.08038	-1.08038	-1.08038	-1.08038	-1.08038	-19556.628059	-19556.628059	-19556.628059	39077.755755	39077.755755	39077.755755	7696.824606
1	2013-01-02	147185.402151	126067.702352	199146.645036	147185.402151	147185.402151	0.00000	0.00000	0.00000	12720.829993	12720.829993	12720.829993	0.00000	0.00000	0.00000	0.00000	0.00000	0.00000	-20026.196678	-20026.196678	-20026.196678	32747.026671	32747.026671	32747.026671	159906.232143
2	2013-01-03	147265.764616	99946.743600	169893.098270	147265.764616	147265.764616	0.00000	0.00000	0.00000	-12982.775088	-12982.775088	-12982.775088	0.00000	0.00000	0.00000	0.00000	0.00000	0.00000	-39372.332251	-39372.332251	-39372.332251	26389.557164	26389.557164	26389.557164	134282.989528
3	2013-01-04	147346.127081	104366.633624	177803.442323	147346.127081	147346.127081	0.00000	0.00000	0.00000	-4071.150474	-4071.150474	-4071.150474	0.00000	0.00000	0.00000	0.00000	0.00000	0.00000	-24171.406596	-24171.406596	-24171.406596	20100.256122	20100.256122	20100.256122	143274.976607
4	2013-01-05	147426.489546	167946.333811	237009.097336	147426.489546	147426.489546	0.00000	0.00000	0.00000	54803.330328	54803.330328	54803.330328	0.00000	0.00000	0.00000	0.00000	0.00000	0.00000	40832.337596	40832.337596	40832.337596	13970.992732	13970.992732	13970.992732	202229.819874

	ds	trend	yhat_lower	yhat_upper	trend_lower	trend_upper	Año Nuevo	Año Nuevo_lower	Año Nuevo_upper	Año Nuevo (observado)	Año Nuevo (observado)_lower	Año Nuevo (observado)_upper	Batalla de Pichincha	Batalla de Pichincha_lower	Batalla de Pichincha_upper	Batalla de Pichincha (observado)	Batalla de Pichincha (observado)_lower	Batalla de Pichincha (observado)_upper	Carnaval	Carnaval_lower	Carnaval_upper	Día de Navidad	Día de Navidad_lower	Día de Navidad_upper	Día de los Difuntos	Día de los Difuntos_lower	Día de los Difuntos_upper	Día del Trabajo	Día del Trabajo_lower	Día del Trabajo_upper	Independencia de Cuenca	Independencia de Cuenca_lower	Independencia de Cuenca_upper	Independencia de Guayaquil	Independencia de Guayaquil_lower	Independencia de Guayaquil_upper	Primer Grito de Independencia	Primer Grito de Independencia_lower	Primer Grito de Independencia_upper	Primer Grito de Independencia (observado)	Primer Grito de Independencia (observado)_lower	Primer Grito de Independencia (observado)_upper	Viernes Santo	Viernes Santo_lower	Viernes Santo_upper	additive_terms	additive_terms_lower	additive_terms_upper	holidays	holidays_lower	holidays_upper	weekly	weekly_lower	weekly_upper	yearly	yearly_lower	yearly_upper	multiplicative_terms	multiplicative_terms_lower	multiplicative_terms_upper	yhat	y	Error	Abs_Error
0	2017-01-01	244017.815789	113959.343408	185581.702752	244017.815789	244017.815789	-201583.851542	-201583.851542	-201583.851542	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	-94793.483742	-94793.483742	-94793.483742	-201583.851542	-201583.851542	-201583.851542	63817.988123	63817.988123	63817.988123	42972.379677	42972.379677	42972.379677	0.0	0.0	0.0	4459.156645	1998.000	147226.332047	147226.332047
1	2017-01-02	244059.992244	243002.909392	313700.574955	244059.992244	244059.992244	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	35137.901325	35137.901325	35137.901325	0.000000	0.000000	0.000000	-1333.837454	-1333.837454	-1333.837454	36471.738779	36471.738779	36471.738779	0.0	0.0	0.0	279197.893568	415182.950	-135985.056432	135985.056432
2	2017-01-03	244102.168699	213830.003854	288913.689794	244102.168699	244102.168699	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	10532.737457	10532.737457	10532.737457	0.000000	0.000000	0.000000	-19373.952292	-19373.952292	-19373.952292	29906.689748	29906.689748	29906.689748	0.0	0.0	0.0	254634.906155	321085.262	-66450.355845	66450.355845
3	2017-01-04	244144.345153	214050.775152	283570.119461	244144.345153	244144.345153	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	3027.895654	3027.895654	3027.895654	0.000000	0.000000	0.000000	-20348.514956	-20348.514956	-20348.514956	23376.410610	23376.410610	23376.410610	0.0	0.0	0.0	247172.240807	274170.595	-26998.354193	26998.354193
4	2017-01-05	244186.521608	189328.165668	257211.903459	244186.521608	244186.521608	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	-22697.301626	-22697.301626	-22697.301626	0.000000	0.000000	0.000000	-39674.946275	-39674.946275	-39674.946275	16977.644649	16977.644649	16977.644649	0.0	0.0	0.0	221489.219982	226863.181	-5373.961018	5373.961018
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
222	2017-08-11	253380.988772	172973.393373	246810.935212	250824.915984	255845.533049	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	-10910.692867	-10910.692867	-10910.692867	0.0	0.0	0.0	0.0	0.0	0.0	-44951.492846	-44951.492846	-44951.492846	-10910.692867	-10910.692867	-10910.692867	-24027.011310	-24027.011310	-24027.011310	-10013.788669	-10013.788669	-10013.788669	0.0	0.0	0.0	208429.495926	238221.922	-29792.426074	29792.426074
223	2017-08-12	253423.165226	247985.978353	317783.916980	250849.791581	255883.908341	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	-970.225589	-970.225589	-970.225589	0.0	0.0	0.0	0.0	0.0	0.0	29886.394697	29886.394697	29886.394697	-970.225589	-970.225589	-970.225589	40940.274164	40940.274164	40940.274164	-10083.653879	-10083.653879	-10083.653879	0.0	0.0	0.0	283309.559923	238627.633	44681.926923	44681.926923
224	2017-08-13	253465.341681	267532.386834	342311.000077	250866.620889	255925.908867	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	53715.995022	53715.995022	53715.995022	0.000000	0.000000	0.000000	63817.988123	63817.988123	63817.988123	-10101.993101	-10101.993101	-10101.993101	0.0	0.0	0.0	307181.336703	258442.452	48738.884703	48738.884703
225	2017-08-14	253507.518136	207551.019493	278110.124063	250883.450198	255998.530756	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	-11394.411615	-11394.411615	-11394.411615	0.000000	0.000000	0.000000	-1333.837454	-1333.837454	-1333.837454	-10060.574161	-10060.574161	-10060.574161	0.0	0.0	0.0	242113.106521	225355.694	16757.412521	16757.412521
226	2017-08-15	253549.694591	188441.118839	258019.785999	250906.585067	256071.152645	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.0	0.000000	0.000000	0.000000	0.0	0.0	0.0	0.0	0.0	0.0	-29327.134620	-29327.134620	-29327.134620	0.000000	0.000000	0.000000	-19373.952292	-19373.952292	-19373.952292	-9953.182328	-9953.182328	-9953.182328	0.0	0.0	0.0	224222.559971	224208.125	14.434971	14.434971

	ds	trend	yhat_lower	yhat_upper	trend_lower	trend_upper	Año Nuevo	Año Nuevo_lower	Año Nuevo_upper	additive_terms	additive_terms_lower	additive_terms_upper	dcoilwtico	dcoilwtico_lower	dcoilwtico_upper	extra_regressors_additive	extra_regressors_additive_lower	extra_regressors_additive_upper	extra_regressors_multiplicative	extra_regressors_multiplicative_lower	extra_regressors_multiplicative_upper	first_weekend_of_month	first_weekend_of_month_lower	first_weekend_of_month_upper	holidays	holidays_lower	holidays_upper	multiplicative_terms	multiplicative_terms_lower	multiplicative_terms_upper	weekly	weekly_lower	weekly_upper	yearly	yearly_lower	yearly_upper	yhat	y
0	2017-01-01	243009.442437	158117.714918	218130.077747	243009.442437	243009.442437	-199407.848469	-199407.848469	-199407.848469	-102625.803428	-102625.803428	-102625.803428	-2982.654978	-2982.654978	-2982.654978	-2982.654978	-2982.654978	-2982.654978	0.196008	0.196008	0.196008	0.196008	0.196008	0.196008	-199407.848469	-199407.848469	-199407.848469	0.196008	0.196008	0.196008	56649.880663	56649.880663	56649.880663	43114.819356	43114.819356	43114.819356	188015.433801	1998.000
1	2017-01-02	243057.520536	249112.745986	305267.447370	243057.520536	243057.520536	0.000000	0.000000	0.000000	34671.418629	34671.418629	34671.418629	-3022.180414	-3022.180414	-3022.180414	-3022.180414	-3022.180414	-3022.180414	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	812.855260	812.855260	812.855260	36880.743782	36880.743782	36880.743782	277728.939165	415182.950
2	2017-01-03	243105.598635	225603.665123	282973.100449	243105.598635	243105.598635	0.000000	0.000000	0.000000	10832.294350	10832.294350	10832.294350	-3140.756722	-3140.756722	-3140.756722	-3140.756722	-3140.756722	-3140.756722	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	-16649.361455	-16649.361455	-16649.361455	30622.412527	30622.412527	30622.412527	253937.892985	321085.262
3	2017-01-04	243153.676734	215599.695640	275729.050398	243153.676734	243153.676734	0.000000	0.000000	0.000000	4199.847249	4199.847249	4199.847249	-2987.204668	-2987.204668	-2987.204668	-2987.204668	-2987.204668	-2987.204668	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	-17244.640661	-17244.640661	-17244.640661	24431.692578	24431.692578	24431.692578	247353.523983	274170.595
4	2017-01-05	243201.754834	194256.315162	249681.638822	243201.754834	243201.754834	0.000000	0.000000	0.000000	-21768.031111	-21768.031111	-21768.031111	-2900.191837	-2900.191837	-2900.191837	-2900.191837	-2900.191837	-2900.191837	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	-37264.848722	-37264.848722	-37264.848722	18397.009448	18397.009448	18397.009448	221433.723723	226863.181

	ds	trend	yhat_lower	yhat_upper	trend_lower	trend_upper	Primer Grito de Independencia	Primer Grito de Independencia_lower	Primer Grito de Independencia_upper	additive_terms	additive_terms_lower	additive_terms_upper	dcoilwtico	dcoilwtico_lower	dcoilwtico_upper	extra_regressors_additive	extra_regressors_additive_lower	extra_regressors_additive_upper	holidays	holidays_lower	holidays_upper	weekly	weekly_lower	weekly_upper	yearly	yearly_lower	yearly_upper	yhat	y
1679	2017-08-11	253682.780473	183103.213852	240088.119375	250375.145048	257195.907307	-6151.014208	-6151.014208	-6151.014208	-42246.608651	-42246.608651	-42246.608651	-3746.434270	-3746.434270	-3746.434270	-3746.434270	-3746.434270	-3746.434270	-6151.014208	-6151.014208	-6151.014208	-20755.577477	-20755.577477	-20755.577477	-11593.582696	-11593.582696	-11593.582696	211436.171822	238221.922
1680	2017-08-12	253730.858573	247203.618742	303393.766577	250401.186149	257256.666668	869.077314	869.077314	869.077314	19960.462377	19960.462377	19960.462377	-3815.817050	-3815.817050	-3815.817050	-3815.817050	-3815.817050	-3815.817050	869.077314	869.077314	869.077314	34451.692392	34451.692392	34451.692392	-11544.490278	-11544.490278	-11544.490278	273691.320950	238627.633
1681	2017-08-13	253778.936672	267281.146361	322439.158588	250434.034626	257316.016968	0.000000	0.000000	0.000000	41353.979025	41353.979025	41353.979025	-3885.199831	-3885.199831	-3885.199831	-3885.199831	-3885.199831	-3885.199831	0.000000	0.000000	0.000000	56649.880663	56649.880663	56649.880663	-11410.701807	-11410.701807	-11410.701807	295132.915697	258442.452
1682	2017-08-14	253827.014771	208299.996618	269736.108706	250456.212196	257396.028471	0.000000	0.000000	0.000000	-14330.313485	-14330.313485	-14330.313485	-3954.582611	-3954.582611	-3954.582611	-3954.582611	-3954.582611	-3954.582611	0.000000	0.000000	0.000000	812.855260	812.855260	812.855260	-11188.586134	-11188.586134	-11188.586134	239496.701287	225355.694
1683	2017-08-15	253875.092870	193961.741771	251148.606134	250478.389766	257461.342467	0.000000	0.000000	0.000000	-31484.719504	-31484.719504	-31484.719504	-3957.994879	-3957.994879	-3957.994879	-3957.994879	-3957.994879	-3957.994879	0.000000	0.000000	0.000000	-16649.361455	-16649.361455	-16649.361455	-10877.363170	-10877.363170	-10877.363170	222390.373367	224208.125

	ds	y	yhat	Error	Abs_Error
1449	2016-12-23	425386.443	302332.025577	-123054.417423	123054.417423
724	2014-12-28	193133.453	312079.894779	118946.441779	118946.441779
1448	2016-12-22	400889.065	287291.374944	-113597.690056	113597.690056
1243	2016-05-31	322106.164	213936.052925	-108170.111075	108170.111075
1214	2016-05-02	388655.729	280896.691823	-107759.037177	107759.037177
720	2014-12-23	357665.078	254897.898205	-102767.179795	102767.179795
361	2013-12-29	186041.000	286952.776135	100911.776135	100911.776135
878	2015-05-31	345804.682	264312.639428	-81492.042572	81492.042572
1084	2015-12-23	367851.184	286855.069760	-80996.114240	80996.114240
1447	2016-12-21	387603.442	307951.005366	-79652.436634	79652.436634
669	2014-11-02	203896.842	282392.555447	78495.713447	78495.713447
1152	2016-03-01	281908.060	204800.254039	-77107.805961	77107.805961
1427	2016-12-01	289689.792	212875.319421	-76814.472579	76814.472579
971	2015-09-01	270499.490	193942.840766	-76556.649234	76556.649234
1446	2016-12-20	384623.784	308666.338852	-75957.445148	75957.445148
1205	2016-04-23	326432.713	250738.384141	-75694.328859	75694.328859
356	2013-12-23	326975.000	252863.726949	-74111.273051	74111.273051
1083	2015-12-22	362154.403	288258.997700	-73895.405300	73895.405300
1001	2015-10-01	248701.812	176229.052031	-72472.759969	72472.759969
1398	2016-11-02	294122.549	223551.109072	-70571.439928	70571.439928
360	2013-12-28	200073.000	269883.579093	69810.579093	69810.579093
1185	2016-04-03	392994.414	324763.936125	-68230.477875	68230.477875
1003	2015-10-03	360545.262	292996.905476	-67548.356524	67548.356524
1305	2016-08-01	287706.632	221547.677507	-66158.954493	66158.954493
1088	2015-12-28	223897.335	289729.682760	65832.347760	65832.347760
1183	2016-04-01	264279.916	200330.706608	-63949.209392	63949.209392
1366	2016-10-01	376487.399	314323.631882	-62163.767118	62163.767118
909	2015-07-01	245410.095	183443.449470	-61966.645530	61966.645530
1212	2016-04-30	334747.624	272907.532606	-61840.091394	61840.091394
1004	2015-10-04	377564.834	315903.104662	-61661.729338	61661.729338

	ds	y	first_weekend_of_month	dcoilwtico	onpromotion	log_onpromotion
0	2013-01-01	700.0	0	0.845816	-0.620444	-1.152637
1	2013-01-02	202020.0	0	0.845816	-0.620444	-1.152637
2	2013-01-03	144878.0	0	0.839389	-0.620444	-1.152637
3	2013-01-04	135754.0	0	0.845060	-0.620444	-1.152637
4	2013-01-05	188356.0	1	0.846068	-0.620444	-1.152637

	ds	y	first_weekend_of_month	dcoilwtico	onpromotion	log_onpromotion
1457	2017-01-01	1998.000	1	-0.660957	-0.574321	-0.640921
1458	2017-01-02	415182.950	0	-0.669715	2.985932	1.841588
1459	2017-01-03	321085.262	0	-0.695992	2.925167	1.830465
1460	2017-01-04	274170.595	0	-0.661965	3.960371	1.998351
1461	2017-01-05	226863.181	0	-0.642683	1.264740	1.418798

	changepoint_prior_scale	seasonality_prior_scale	holidays_prior_scale	seasonality_mode	n_changepoints	yearly_seasonality	weekly_seasonality	rmse
2639	0.010	0.1	0.05	multiplicative	10	30	7	41878.212478
1391	0.005	5.0	0.05	multiplicative	50	30	7	41896.214139
3295	0.010	1.0	0.05	multiplicative	25	30	7	41899.175549
2655	0.010	0.1	0.05	multiplicative	25	30	7	41904.746282
2031	0.005	10.0	0.05	multiplicative	50	30	7	41905.102645
1999	0.005	10.0	0.05	multiplicative	10	30	7	41905.790007
4559	0.010	10.0	0.05	multiplicative	10	30	7	41908.741257
3279	0.010	1.0	0.05	multiplicative	10	30	7	41909.878317
1519	0.005	5.0	0.10	multiplicative	50	30	7	41911.982345
4607	0.010	10.0	0.05	multiplicative	75	30	7	41912.646285
719	0.005	1.0	0.05	multiplicative	10	30	7	41913.349810
2015	0.005	10.0	0.05	multiplicative	25	30	7	41914.081858
895	0.005	1.0	0.10	multiplicative	75	30	7	41916.502421
2175	0.005	10.0	0.10	multiplicative	75	30	7	41920.932640
207	0.005	0.1	0.10	multiplicative	10	30	7	41921.397098
3311	0.010	1.0	0.05	multiplicative	50	30	7	41923.662296
1375	0.005	5.0	0.05	multiplicative	25	30	7	41926.530679
3951	0.010	5.0	0.05	multiplicative	50	30	7	41927.321722
2671	0.010	0.1	0.05	multiplicative	50	30	7	41927.974469
735	0.005	1.0	0.05	multiplicative	25	30	7	41928.895116

	horizon	mse	rmse	mae	mape	mdape	smape	coverage
0	23 days	2.804676e+09	52959.195092	25410.413904	0.088853	0.054090	0.097145	0.684956
1	24 days	1.942751e+09	44076.648208	22749.639352	0.083952	0.054163	0.090155	0.693805
2	25 days	1.570253e+09	39626.423221	21020.516048	0.080114	0.054163	0.084912	0.697345
3	26 days	1.169486e+09	34197.751055	19039.718365	0.074903	0.053285	0.078149	0.715044
4	27 days	7.522528e+08	27427.226708	17293.462739	0.070621	0.052459	0.072122	0.720354
...	...	...	...	...	...	...	...	...
200	223 days	8.956983e+08	29928.218533	23520.973729	0.099555	0.092779	0.093899	0.437168
201	224 days	9.234287e+08	30387.969933	23905.742922	0.101677	0.091499	0.096409	0.431858
202	225 days	9.207709e+08	30344.206579	23924.946763	0.102949	0.087830	0.097492	0.438938
203	226 days	1.016762e+09	31886.704821	25320.287895	0.108900	0.092779	0.102783	0.419469
204	227 days	1.086775e+09	32966.272665	26405.895833	0.114392	0.094658	0.107456	0.394690

Time-series Forecasting Using Facebook Prophet¶

Alan Gaugler¶

December 2024¶

1. Introduction¶

2. Description of the Dataset¶

Main Files¶

Supplementary Files¶

Additional Notes¶

Key Observations¶

Train-Test Split¶

3. Project setup¶

Import the Required Libraries¶

Create Various Functions¶

Load the datasets¶

4. Exploratory Data Analysis and Data Preparation¶

Plotting and Aggregating Sales Data by Product Category¶

Plot characteristics of the dataset¶

Data Cleaning¶

Outlier Treatment¶

5. Modelling¶

Prepare the dataset for modelling with Prophet¶

Prophet Model 1. Basic Prophet Model¶

Determine the number of days in the test set¶

Plot and Analyze the Forecast¶

Plot Observations:¶

Observations from the Plot:¶

Calculate the Performance Metrics for the Training Set: Model Pr1¶

Performance Metrics:¶

Insights:¶

Plot the Test Set Predictions in More Detail¶

Observations¶

Calculate the Performance Metrics on the Test Set: Model Pr1¶

Display the 5 days with the largest prediction errors¶

Key Insights¶

Prophet Model 2.1. Incorporating Holidays¶

Observations: Changes in the Forecasted DataFrame After Introducing Holiday Regressors¶

Key Differences After Adding Holiday Regressors:¶

Observations from the Plot¶

Breakdown the components of the series¶

Plot and Analyze the Forecast¶

Plot Predicted vs. Actual Sales¶

Observations from the Plot¶

Calculate Performance Metrics for the Training Set: Model Pr2.1¶

Plot the test set predictions in more detail¶

Calculate Performance Metrics for the Test Set: Model Pr2.1¶

Prophet Model 2.2. Model with New Year's As A Multiplicative Holiday¶

Create Flags for Multiplicative Holidays on the Main Grocery Dataframe¶

Prepare the data frame for modelling with Prophet¶

Plot and Analyze the Forecast¶

Plot Predicted vs. Actual Sales¶

Calculate Model Performance Metrics for the Training Set: Model Pr2.2¶

Calculate Model Performance Metrics for the Test Set: Model Pr2.2¶

Prophet Model 2.3. Create a prophet model for New Year's Day Only and Merge with Pr2.1¶

Options Considered¶

Chosen Solution¶

Create a Function for New Year's Day Only and Merge with Pr2.1¶

Plot and Analyze the Forecast¶

Calculate Performance Metrics for the Test Set: Model Pr2.3¶

Display the 10 Days with the Largest Prediction Errors¶

Key Insights from Model 2.3¶

Model 3.1. Additional Regressors Including Oil Prices, First Weekend After Payday and Earthquake Outliers¶

Existing Regressors¶

Observations from Model Pr2.1¶

New Regressors¶

Extreme Outliers¶

Plot Daily Sales for April 2016¶

Define and Add Key Earthquake Dates as Regressors¶

Log Transformation of the onpromotion Variable¶

Prepare the data frame for modelling with Prophet¶

Model Configuration with Holiday Data and Additional Regressors¶

Select the Additional Regressors¶

Plot and Analyze the Forecast¶

Plot Predicted vs. Actual Sales¶

Observations from the Plot¶

Comparison to Model Pr2.2¶

Calculate Performance Metrics for the Training Set: Model Pr3.1¶

Create a Model for New Year's Day Only and Merge with Pr3.1¶

Observations from the Plot¶

Comparison to Model Pr2.3¶

Calculate Performance Metrics for the Test Set: Model Pr3.1¶

Log Transformation of the `onpromotion` Variable¶

Model	RMSE	MAE	MAPE	Accuracy (%)	Key Additions
Pr1	41,012.7	24,738.0	81.51%	18.49%	Baseline model with general seasonality and trend components.
Pr2.1	33,199.8	21,876.6	40.20%	59.80%	Incorporated holiday effects into the model.
Pr2.2	32,534.50	21,782.7	25.66%	74.34%	Treated New Year’s as a multiplicative holiday. Not implemented as Pr2.3's solution was far better.
Pr2.3	31,729.6	21,238.8	8.28%	91.72%	Created a separate model for New Year’s Day.
Pr3.1	28,133.3	19,643.9	7.76%	92.24%	Added regressors for oil prices, weekends after payday and the outlier impact of the earthquake.
Pr3.2	26,952.3	18,853.5	7.49%	92.51%	Added regressors for the 1st day of the month and 28 December.
Pr4	26,955.3	18,852.6	7.49%	92.51%	Scaled non-binary regressors for consistent influence. This had no effect so was not implemented.
Pr5	24,679.2	16,990.6	7.03%	92.97%	Applied cross-validation and hyperparameter tuning.